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SUMMARY 


Developments,  applications,  and  evaluation  of  the  SCIPUFF  (Second-order 
Closure  Integrated  PUFF)  model  are  presented.  SCIPUFF  is  a  Lagrangian  puff  model  for 
atmospheric  transport  and  diffusion  of  nuclear  clouds  at  late-time,  that  is,  times  of  hours 
after  initiation.  The  model  uses  turbulence  closure  theory  to  relate  the  atmospheric 
dispersion  of  the  cloud  to  measurable  statistics  of  the  velocity  field,  and  also  contains  a 
prediction  of  the  'uncertainty'  in  a  model  prediction  arising  from  the  random  velocity 
fluctuations. 

SCIPUFF  has  been  extended  to  iiK:lude  a  complete  tensor  description  of  the 
second-order  spatial  moments,  allowing  an  accurate  representation  of  the  shear-induced 
distortion  of  a  cloud.  The  model  also  incorporates  the  observational  spectra  of  Nastrom 
and  Gage  (1985)  for  the  large  scale  velocity  fluctuations,  and  uses  this  information  to 
describe  the  diffusion  of  a  cloud  and  also  to  characterize  the  uncertainty  in  a  finite- 
resolution  representation  of  the  large  scale  wind  field.  An  option  to  use  a  'relative' 
diffusion  prediction  has  been  implemented  to  facilitate  comparison  with  deterministic 
dispersion  models.  This  approximation  neglects  the  meandering  effect  of  eddies  on 
scales  larger  than  the  cloud,  which  are  effectively  assumed  to  contribute  only  to  an 
uncertainty  in  position. 

An  extensive  evaluation  exercise  has  been  conducted  using  the  Across  North 
America  Tracer  Experiment  (ANATEX).  This  experiment  collected  a  large  data  base  of 
high  quality  tracer  concentration  measurements  on  continental  scales  for  a  three  month 
period.  Surface  dose  samples  and  short  term  aircraft  samples  were  predicted  and 
reasonably  good  agreement  obtained.  The  wind  fields  are  not  accurate  enough  to  allow  a 
reliable  prediction  of  point  concentrations,  but  the  statistical  distribution  agreed  well  with 
observations.  The  aircraft  data,  in  particular,  required  the  correct  estimate  of  the 
probability  distribution  to  match  the  occasional  high  concentrations  in  the  instantaneous  — 
plume.  Some  skill  was  also  demonstrated  in  a  prediction  based  on  a  'climatological' 
average  of  the  wind  field,  .i.e.,  no  detailed  information  about  the  flow. 

;d  □ 

A  number  of  numerical  integrations  have  been  made  to  examine  the  effects  of  . : 

atmospheric  boundary  layer  turbulence  on  the  evolution  of  a  nuclear  dust  cloud  over . . 

periods  of  the  order  of  a  day.  A  particular  concern  is  the  fate  of  the  low-level  dust  ^ _ 

pedestal,  containing  high  concentrations  of  material  swept  up  from  the  surface,  which  Codes 


j  rtvdii  and  /  or 


may  be  lofted  by  turbulent  eddies  and  remain  in  the  atmosphere  for  long  periods. 
Calculations  have  been  made  for  a  range  of  atmospheric  conditions  for  periods  of  24 
hours.  The  general  conclusion  from  these  studies  is  that  the  major  factor  influencing  die 
fraction  of  dust  remaining  aloft  after  12-24  hours  is  the  initial  size  distribution. 
Meteorological  variations  produce  some  effect  on  the  rate  of  fallout,  but  material  with 
diameters  above  about  30pm  is  mostly  deposited  within  12  hours,  while  diameters  below 
10pm  are  only  slightly  depleted  by  gravitational  settling  over  this  time  period.  The 
turbulent  deposition  can  be  significant  over  forested  areas,  however,  so  that  particles  of 
10pm  size  can  be  significantly  reduced  over  a  12  hour  period  when  compared  to  1pm 
particles. 

An  analytical  description  of  the  statistics  of  low-level  dust  concentrations  at  late 
time  after  a  nuclear  burst  was  obtained,  based  on  the  r.m.s.  fluctuations  in  the  wind  field. 
The  problem  is  complicated  by  the  fact  that  late-time  dust  at  low  altitude  is  comprised  of 
a  superposition  of  different  particle  sizes  falling  from  different  initial  altitudes.  The 
spatial  location  of  the  different  size  groups  therefore  depends  on  the  integrated  wind 
vector  over  the  range  of  the  descent  The  analytic  predictions  were  compared  with  a 
direct  Monte  Carlo  simulation  of  a  horizontally  homogeneous  wind  field  acting  on  the 
initial  cloud.  The  analytic  description  must  account  for  the  vertical  correlation  of  the 
wind  field.  A  simple  exponential  correlation  function  resulted  in  a  good  representation  of 
the  probability  distribution  for  an  aircraft  intercept  mass. 

A  study  was  performed  to  determine  the  sensitivity  of  the  dispersion  prediction  to 
variations  in  the  initial  conditions  for  a  release  of  radioactive  material  from  an  attack  on  a 
nuclear  power  plant.  The  diffusive  nature  of  atmospheric  dispersion  reduces  the 
dependence  on  initial  conditions  as  transport  time  increases,  although  conserved 
quantities  such  as  total  release  mass  are  clearly  unaffected  by  the  diffusion  process. 
Factors  such  as  release  height  and  release  size  are  not  critical  for  transport  over  distances 
greater  than  a  few  kilometers.  The  most  critical  aspects  ate  the  release  height  relative  to 
the  planetary  boundary  layer  depth  and  the  size  distribution  of  the  release.  Releases  of 
small  particles  or  gaseous  species  above  the  mixed  layer  can  remain  aloft  for  long  times 
before  being  mixed  down  to  the  surface.  Releases  inside  the  mixed  layer,  however,  are 
quickly  mixed  through  the  boundary  layer.  A  demonstration  calculation  was  performed 
using  a  DICE/MAZ  prediction  of  the  radioactive  cloud  from  a  20kT  we:q>on,  showing  the 
extmit  of  the  contamination  and  the  probabilistic  prediction  capability  of  the  model. 
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CONVERSION  TABLE 

Conversion  Factors  for  U.S.  Customary  to  Metric  (SI)  Units  of  Measure 


MULTIPLY 

TO  GET 

BY 

BY 

TO  GET 

DIVIDE 

angstrom 

1.000  000  X  E -10 

meters  (m) 

atnwsphere  (normal) 

1.013  25  XE-t-2 

kilo  pascal  (kPa) 

bar 

1.000  000  XE  +2 

kiio  pascal  (kPa) 

British  thermal  unit  (thermomechanical} 

1.054  350  XE  43 

joule  ( J) 

calorie  (thermomechanicai) 

4.184  000 

joule  (J) 

cal  (thermochemical)/cm^ 

4.184  000  000  X  E -2 

megajoula/m^  (MJ/M^) 

curie 

3.700  000  X  4- 

*  giga  becquerei  (GBq) 

degree  (angle) 

1.745  329  XE -2 

radian  (rad) 

degree  Fahrenheit 

tk* **(TF4-459.67)/1.8 

degree  kelvin  (K) 

electron  volt 

1.602  19  XE -19 

joule  ( J) 

erg 

joule  (J) 

erg^second 

1.000  000  XE -7 

watt  (W) 

foot 

3.048  000  X  E  -1 

meter  (m) 

foot-pound-force 

1.355  818 

joulie  (J) 

gallon  (US.  liquid) 

3.785  412  X  E  -3 

meter^  (nP) 

inch 

2.540  OOOXE-2 

meter  (m) 

jerk 

1.000  000  XE  4-9 

joule  (J) 

joule/tdlogram  (J/kg)  (radiation 
dose  absorbed) 

1.000  000 

Gray  (Gy) 

kilotons 

4.183 

terajouies 

kip  (1000  fcf) 

4.448  222XE43 

newton  (N) 

kip/inch^  (ksi) 

6.894  757  XE  43 

kilo  pascal  (kPa) 

ktap 

1.000  000  X  E  4-2 

newton-seconcym^ 

(N-a/m2) 

micron 

1.000  OOOXE-e 

meter  (m) 

mil 

2.540  000  X  E  -5 

meter  (m) 

mile  (international) 

1.609  344  XE  43 

meter  (m) 

ounce 

2.834  952XE-2 

kilogram  (kg) 

pound-force  (lbs  avoirdupois) 

4.448  222 

newton  (N) 

pound-force  inch 

1.129  848  XE'1 

newton-meter  (N*m) 

pound/force/Rich 

1.751  288  XE  4-2 

newton/meter  (N/m) 

poundfforce/foot^ 

4.788  026  X  E  -2 

k9o  pascal  (kPa) 

pound-force  inch^  (p*0 

6.894757 

kiio  pascal  (kPa) 

pound-mass  (bm  avoirdupois) 

4.535  924  X  E  -1 

kilogram  (kg) 

pound-mass-foot^  (moment  of  inertia) 

4.214  011  XE-2 

kHogranvmeter 

pound-mass-foot^ 

1.601  846  XE  4-1 

kilogranYmeter^  (kg^m?) 

rad  (radiation  dose  absorbed) 

1.000  000  XE-2 

-  Gray  (Gy) 

roentgen 

2.579  760  XE -4 

coulomb/kilogramm  (C/kg) 

shake 

1.000  000  XE -8 

second  (s) 

slug 

1.459  390  X  E  4-1 

kilogrem  (kg) 

torr  (mm  Hg,  0*0) 

1.333  22XE-1 

kilo  pascal  (kPa) 

•  the  becquerel  (Bq)  is  the  SI  unit  of  radioaetiviiy;  1  Bq  =  1  event^s. 

**  The  Grey  (^)  is  the  SI  unit  of  absorbed  recSation. 


V 


TABLE  OF  CONTENTS 


Sectkm  Page 

SUMMARY .  ui 

CONVERSION  TABLE .  v 

HGURES .  vii 

1  INTRODUCTION . .  1 

2  DUST-OFF  MODEL  COMPARISON . .  3 

2.1  OBJECTIVES  AND  METHODOLOGY .  3 

2.2  SCIPUFF  MODEL  IMPROVEMENTS _  3 

2.2.1  Improved  Wind  Shear  Treatment .  4 

2.2.2  Relative  Diffusion  Prediction .  7 

2.3  TEST  CASE  RESULTS . .  9 

3  ANATEX  MODEL  EVALUATION  STUDY . 18 

4  BOUNDARY  LAYER  DIFFUSION  EFFECTS .  27 

4.1  BACKGROUND _ 27 

4.2  NUMERICAL  MODEL . 28 

4.2.1  Turbulent  Transport  Model .  28 

4.2.2  Turbulent  Dqxrsition  Model _  30 

4.2.3  Initial  Conditions .  32 

4.2.4  Meteorological  Parameters . . .  35 

4.3  RESULTS _ 38 

4.4  SUMMARY  OF  BOUNDARY  LAYER  RESULTS .  52 

5  STATISTICAL  FALL-OUT  ANALYSIS _ 53 

6  SCIPUFF  MODEL  STUDIES _ 69 

6.1  RADIOACTIVE  FALLOUT  FROM  20KT  WEAPON _  69 

6.2  NUCLEAR  POWER  PLANT  RELEASE  SENSITTVmES .  77 

7  CONCLUDING  REMARKS . 86 

8  REFERENCES _  89 

APPENDIX . 91 

vi 


FIGURES 


Figure  Page 

2-1  Multibuist  scenario  release  locations .  4 

2-2  Winds  for  16  January  1978  at  SS^N,  37.5<»E  and  OOZ .  5 

2-3  SCIPUFF  results  for  the  single  burst  scenario  using  the  January  wind 
field.  Perspective  views  of  the  10'*°  g/cm^  iso-surface  at  r  =:  0,  1  and  8 
hrs .  8 

2-4  SCIPUFF  results  for  the  single  burst  scenario  using  the  January  wind 
field  Horizontal  cross-sections  at  r  =:  8  hrs.  Contour  levels  of  10~^^ 

10-12. 10-11, 10-10, 10-9,  i(H,  and  lO-^  gm/cm^ .  11 

2-S  SCIPUFF  results  for  the  single  burst  scenario  using  the  July  wind  field. 
Horizontal  cross-sections  at  r  =  8  hrs.  Contour  levels  of  lO-i^,  10-12, 

10-11, 10-10, 10-9,  io-«,  and  io-7  gm/cm^ .  12 

2-6  SCIPUFF  results  for  the  multi-burst  scenario  using  the  January  wind 
field.  Horizontal  cross-sections  at  r  ^  g  hrs.  Contour  levels  of  lO-i^, 

10-12, 10-11. 10-10, 10-0, 10-8,  and  10-2  gm/cm^ .  13 

2-7  SCIPUFF  results  for  the  multi-burst  scenario  using  the  July  wind  field. 
Horizontal  cross-sections  at  r  =  8  hrs.  Contour  levels  of  10*1^  10'i2, 10* 

11, 10-10,  lO-o,  10-»,  and  10-2  gm/cm^ .  14 

2-8  Intercepted  dust  masses  along  horizontal  flight  paths  through  the  single 
burst  cloud  after  8  hours  with  winter  winds,  (a)  East-West  paths, 

(b)  North-South  paths.  - z=10(hn; - z=8km; 

- z=12km .  15 

2-9  SCIPUFF  results  for  the  single  burst  scenario  using  the  January  wind 
field  and  absolute  diffusion  scheme.  Horizontal  cross-sections  at  r  =  8  hrs 
and  z  =  100  m.  (a)  mean  total  dust  concentration,  (b)  r.m.s.  fluctuation  in 
total  dust  Contour  levels  of  10-l^  IO-12,  lO-n,  lO-i®,  and  10"®  gm/cm^ .  16 


vu 


FIGURES  (C<Mitinued) 


Figure  Page 

2- 10  Fractal  realizations  from  the  ensemble  statistics  of  Figure  2-9.  Horizontal 

cross-sections  at  r  =  8  hrs  and  z=100m.  Contour  levels  of  lO"^^, 

10“^^  10“^®,  and  10^  gm/cm^ .  16 

3- 1  ANATEX  release  and  surface  sampler  locations .  19 

3-2  Comparison  of  the  SCIPUFF  24hr  avmge  FTCH  surface  concentration 

patterns  with  the  observed  patterns  for  tl^  period  23-25  March .  21 

3-3  Predicted  and  observed  ground  sample  cumulative  distribution  functions 

for  the  period  9-28  March .  22 

3-4  Predicted  and  observed  ground  sample  cumulative  distribution  functions 
for  PTCH  along  three  arcs  for  the  period  5  January-21  February. 
Predictions  were  made  using  a  standard  boundary  layer  depth  of  1000m .  23 

3-S  Predicted  and  observed  aircraft  sample  cumulative  distribution  functions 

for  the  period  9-28  March .  24 

3-6  Adjusted  predicted  and  observed  aircraft  sample  cumulative  distribution 

functions  for  the  period  9-28  March .  25 

3- 7  Predicted  and  observed  ground  sample  cumulative  distribution  functions 

for  the  period  9-28  March.  Predictions  were  made  using  the  ensemble 
meteorology .  26 

4- 1  Surface  deposition  and  gravitational  settling  velocity  versus  particle 

diameter  for  two  surface  roughness  lengths.  Surface  friction  velocity  is 
assumed  to  be  0.4m/s,  and  the  deposition  is  calculated  for  a  height  of 
10m .  32 

4-2  Dust  field  from  DICE  at  r  =  5  mins.  Contours  are  in  g/cc .  33 


HGUPES  (Continued) 


flsure  Page 

4-3  Vertical  profiles  of  horizontally-integrated  dust  concentration  from  the 
DICE  run  at  r  =  5  mins,  (a)  show  full  profile  (note  logarithmic 
concentration  scale),  (b)  shows  detail  in  the  lowest  levels .  34 

4-4  Diurnal  variation  of  the  planetary  boundary  layer,  (a)  r.m.s.  vertical 
velocity  fluctuation,  contour  interval  O.lm/s.  (b)  surface  heat  flux  in 
Km/s .  37 

4-5  Integrated  dust  profiles  at  4,  8,  and  12  hours  after  release.  Release  is  at 
(a)  ()0:(X);  (b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  ^)eed  lOm/s,  roughness 
length  0.01m . . .  39 

4-6  Integrated  dust  profiles  at  4,  8.  and  12  hours  after  release.  Release  is  at 
(a)  00:00;  (b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  lOm/s,  roughness 
length  1.0m .  42 

4-7  Integrated  dust  profiles  at  4,  8,  and  12  hours  after  release.  Release  is  at 
(a)  00:(X);  (b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  2m/s,  roughness 
length  1.0m .  44 

4-8  Integrated  dust  profiles  at  12  hours  after  06:00  release.  Surface  heat  flux 

maximum  is  0.075Km/s  unless  noted .  46 

4-9  Integrated  dust  profiles  at  4,  8,  and  12  hours  after  release.  Release  is  at 
(a)  00:00;  (b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  qieed  lOm/s,  roughness 
length  0.01m,  pedestal  dust  only . 48 

4-10  Integrated  dust  profiles  at  12  hours  after  06:00  release  assuming  all  dust 

is  10  microns.  Surface  heat  flux  maximum  is  0.07SKm/s  unless  noted .  50 

4-1 1  Integrated  dust  profiles  at  12  hours  after  06:00  release  assuming  proposed 
incipient  particle  size  distribution.  Surface  heat  flux  maximum  is 
0.075Km/s  unless  noted . 51 

ix 


FIGURES  (Continued) 


Figure  Page 

S-1  Wind  statistics  for  the  East-West  component  from  3652  Moscow  wind 
profiles,  (a)  local  wind  at  altitude;  (b)  average  wind  between  altitude  and 
ground .  57 

5-2  Mean  mass  loading  after  8  hours  along  a  North-South  flighq)ath  at  100m 
altitude  for  an  ensemble  of  3652  Moscow  winds.  The  explicit  Monte- 
Carlo  result  (solid)  and  the  ensemble  statistical  calculation  (dashed)  ate 
shown  as  a  fimction  of  East-West  location  relative  to  the  burst .  58 

5-3  Cumulative  distribution  of  the  mass  loading  at  8  hours  from  the  Monte- 
Carlo  results  with  3652  winds  (i.e.  probability  that  loading  is  less  than 
given  value)  at  two  locations .  60 

5-4  Correlation  between  East-West  velocities  at  two  altitudes  from  3652 
Moscow  wind  profiles.  Lower  level  is  fixed  at  14S0m.  Correlations 
between  the  local  velocity  and  the  avoage  over  the  layer  between  the  two 
levels  are  shown .  63 

5-5  R.m.s.  mass  loading  after  8  hours  for  a  North-South  fligh^ath  at  100m 
from  an  ensemble  of  3652  Moscow  winds.  The  explicit  Monte-Carlo 
result  (solid)  and  the  statistical  ensemble  calculation  (dashed)  are  shown .  65 

5-6  Probability  of  mass  loading  exceeding  1  mg/cm^  for  the  ensemble  of 
3652  Moscow  winds.  Comparison  between  Monte-Carlo  calculation 
(solid)  and  statistical  estimate  (dashed) _  66 

5- 7  Probability  of  mass  loading  exceeding  10  mg/cm^  for  the  ensemble  of 

3652  Moscow  winds.  Comparison  between  Monte-Carlo  calculation 
(solid)  and  statistical  estimate  (dashed) .  67 

6- 1  Distribution  of  radioactive  material  at  r  =:  5mins.  (a)  vertical  distribution 

of  active  mass;  (b)  distribution  by  particle  dze,  contours  indicate  fraction 
of  mass  contained  in  particles  smaller  than  the  indicated  diameter. .  70 


X 


FIGURES  (Continued) 


F^ure  Page 

6-2  Initial  wind  field  for  the  SCIPUFF  simulation .  71 

6-3  SCIPUFF  prediction  of  the  ensemble  mean  concentration  and  standard 
deviation  at  z=100m  and  z=5km  8  hours  after  initiation  of  a  20kT  burst 
centered  at  42*N,  88*W.  Contour  levels  of  10"2,  10"^,  1,  and 
lOpg/m^ . 73 

6-4  SCIPUFF  prediction  of  the  ensemble  mean  surface  dose  8  hours  after 
initiation  of  a  20kT  burst  centered  at  42*N,  88‘W.  Contour  levels  of  10"^, 

10^,  10"5, 1(H,  1(H  and  10~2)ig/m2 .  74 

6-S  SCIPUFF  prediction  of  the  ensemble  mean  surface  radiation  dose  8  hours 
after  initiation  of  a  20kT  burst  centered  at  42*N,  88*W.  Contour  levels  of 
0.01, 0.1, 10-5, 1, 10  and  lOOrad/hr .  75 

6-6  SCIPUFF  prediction  of  the  probability  of  exceeding  a  surface  radiation 
dose  of  (a)  O.lrad/hr  and  (b)  lOrad/hr  8  hours  after  initiation  of  a  20kT 
burst  Contour  levels  of  1%,  5%,  20%,  50%,  80%,  100% .  76 

6-7  Surface  dose  at  48  hours  after  the  release  of  a  gaseous  material  at  500m 
from  the  Illinois  source.  Release  start  at  (a)  01/07/87:00:00, 

(b)  02/04/87:00:00.  (c)  03/04/87:00:00,  (d)  01/07/87:12:00, 

(e)  02/04/87:12:00,  (003/04/87:12:00 .  79 

6-8  Surface  dose  at  48  hours  after  the  release  of  a  gaseous  material  at  5(X)m 
from  the  Four  Comers  source.  Release  start  at  (a)  01/07/87 :(X):00, 

(b)  02/04/87:00:00,  (c)  03/04/87:00:00.  (d)  01/07/87:12:00, 

(e)  02/04/87: 12:00,  (f)  03/04/87:12:00. .  79 

6-9  Surface  dose  at  48  hours  after  the  release  on  02/04/87  at  12:00  firom  the 
Illinois  source.  Release  of  (a)  gas  at  500m,  (b)  gas  at  1000m,  (c)  gas  at 
2000m,  (d)  lOp  at  50Om,  (e)  lOp  at  1000m,  (f)  lOp  at  2000m .  80 


FIGURES  (Continued) 


Figure  Page 

6-10  Surface  dose  at  48  hours  after  the  release  on  02/04/87  at  12:00  from  the 
Illinois  source.  Release  of  (a)  3(Him  at  SOOm,  (b)  30^m  at  lOOOm, 

(c)  30|im  at  2000m .  80 

6-11  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m 
from  the  Illinois  source.  Release  start  at  (a)  01/07/87:00:00, 

(b)  02/04/87:00:00,  (c)  03/04/87:00:00,  (d)  01/07/87:12:00, 

(e)  02/04/87: 12:00,  (f)  03/1^/87:12:00. . 82 

6-12  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m 
from  the  Illinois  source..  Release  start  at  (a)  01/07/87:00:00, 

(b)  02/04/87:00:00,  (c)  03/04/87:00:00,  (d)  01/07/87:12:00, 

(e)  02/04/87: 12:00,  (f)  03/04/87:12:00. .  82 

6-13  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m 
from  the  Sliaois  source,  (a)  gas  at  20Qm,  (b)  gas  at  1000m,  (c)  gas  at 
2000m,  (d)  10^  at  200m,  (e)  IO4  at  1000m,  (f)  KV  at  2000m .  83 

6-14  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m 
from  the  Illinois  source,  (a)  30^m  ^  200m,  (b)  30(tm  at  1000m,  (c) 

30^m  at  2000m,  (d)  lOO^m  at  200m,  (e)  100)im  at  1000m,  (f)  100^m  at 
2000m .  83 

6- IS  Surface  dose  at  8  hours  after  the  release  on  02/04/87  at  00:00  from  the 
Illinois  source,  (a)  gas  at  200m,  (b)  gas  at  SOOm,  (c)  gas  at  1000m, 

(d)  lOiim  at  200m,  (e)  lO^im  at  SOOm,  (f)  lO^im  at  1000m .  84 

6-16  Surface  dose  at  8  hours  after  the  release  on  02/04/87  at  00:00  from  the 
Illinois  source,  (a)  SO^m  at  200m,  (b)  30iim  at  SOOm,  (c)  30itm  at 
1000m,  (d)  lOO^im  at  200m,  (e)  100|xm  at  SOOm,  (f)  lOO^un  at  1000m .  84 


xii 


SECTION  1 


INTRODUCTION 

The  dust  cloud  lofted  by  a  nuclear  explosion  contains  a  wide  range  of  particle 
sizes,  including  many  fine  particles  that  can  remain  airborne  for  hours  or  days  after  the 
detonation.  The  dust  environments  from  multiple  bursts  can  therefore  present  a  hazard  to 
various  defense  systems  and  there  is  a  need  for  quantitative  prediction  of  the  late-time 
cloud  distributions.  The  late-dme  evolution  of  the  cloud  is  largely  controlled  by  the 
atmospheric  wind  field,  which  is  responsible  for  dispersing  the  small  particles.  Other 
factors,  such  as  gravitational  settling  and  interaction  with  water/ice  clouds,  also  play  a 
role  in  determining  the  removal  rate  of  the  dust  firom  the  atmosphere. 

The  effects  of  wind  on  a  dust  cloud  are  generally  described  as  transport  and 
diffusion,  and  most  mathematical  models  make  a  clear  conceptual  distinction  between  the 
two  processes.  The  diffusion  is  represented  in  an  analogous  way  to  molecular  difiusion 
and  the  observed,  or  forecast,  wind  field  is  used  to  move  the  cloud.  This  distinction 
between  transport  and  diffiision  relies  upon  a  wide  difference  in  length  scales  for  the  two 
processes.  The  diffusion  of  the  cloud  is  actually  brought  about  by  random  turbulent 
eddies,  which  separate  particles  from  one  another  as  they  are  advected  in  the  chaotic  flow 
field.  If  the  turbulence  eddies  only  exist  at  small  scales,  then  we  can  clearly  distinguish 
between  the  transport  winds  (the  large  scale  variations)  and  the  diffusing  winds. 
Unfortunately,  the  spectrum  of  ''^Idies  in  the  atmosphere  includes  variations  on  all  spatial 
scales,  from  millimeters  u  ,lobal  scale,  and  there  is  no  spectral  gap  between  the 
resolvable  transport  winds  ano  the  random  turbulent  eddies.  The  distinction  between 
transport  and  diffusion  is  therefore  arbitrary,  and  is  defined  by  the  resolution  of  the  wind 
field.  AU  "known"  details  of  the  wind  field  are  included  as  tran^rt,  while  "unknown" 
wind  variations  must  be  represented  statistically  as  diffusion.  The  "knowledge"  of  the 
wind  field  is  determined  by  the  observing  network,  or  the  reliability  of  the  numerical 
forecast 

The  preceding  discussion  clearly  implies  that  no  prediction  of  atmospheric 
dispersion  can  be  completely  certain.  The  degree  of  uncertainty  depends  on  the  available 
input  information  for  the  prediction,  but  the  inherent  limitation  on  our  knowledge  of  the 
detailed  state  of  the  atmosphere  generally  leads  to  ngnificant  uncertainty  levels.  The 
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initial  development  of  the  SCIPUFF  (Second-order  Closure  Integrated  PufQ  model 
attempted  to  include  a  quantitative  treatment  of  the  uncertainty  due  to  three-dimensional 
turbulent  eddies  in  die  planetary  boundary  layer  (Sykes  et  al.,  1988),  and  was  extensively 
validated  against  local  di^iersion  measurements  at  ranges  of  up  to  50km.  Under  DNA 
sponsorship,  SCIPUFF  was  extended  to  include  dynamic  effects  of  water  condensation 
and  evaporation  (Lewellen  et  aL,1991)  and  also  for  application  to  larger  scale  dispersion 
problems.  In  this  report,  we  describe  the  continued  development  of  the  quantitative 
description  of  prediction  uncertainty  for  long  range  atmospheric  transport  and  dispersion. 

Several  studies  are  reported  below.  &rtensions  of  SCIPUFF  and  comparisons 
widi  other  models  and  with  tracer  data  from  field  experiment  are  described  in  Sections  2 
and  3.  An  investigation  of  the  effects  of  vertical  diffusion  in  the  planetary  boundary  layer 
is  presented  in  Section  4.  Section  5  deals  with  the  representation  of  the  uncertainty  in  a 
falling  cloud  of  dust  particles,  with  application  to  the  aircraft  environment  at  low  altitude 
at  times  of  several  hours  after  burst  Section  6  describes  several  applications  of  SCIPUFF 
in  the  prediction  of  atmospheric  dispersion,  and  finally  Section  7  summarizes  the 
principal  findings  of  the  report 
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SECTION  2 


’DUST-OFF’  MODEL  COMPARISON 

2.1  OBJECTIVES  AND  METHODOLOGY. 

The  'Dust-Off  model  comparison  exercise  was  initiated  to  examine  the 
differences  between  various  late-time  cloud  prediction  models,  and  to  determine  the 
important  features  of  the  models  with  respect  to  the  prediction  of  nuclear  cloud  evolution. 
Current  mottels  use  a  variety  of  diffusion  parameterizations  and  numerical  representation 
techniques,  and  valuable  information  on  sensitivities  can  be  obtained  from  direct 
comparison  of  model  predictions.  A  series  of  test  cases  and  a  set  of  output  measures 
were  defined  for  three  models,  SCIPUFF,  SAICs  TORAS,  and  Kaman's  PEPPER  code. 

The  test  cases  consisted  of  both  single  and  multiburst  scenarios  and  two 
representative  wind  fields.  The  initial  cIoikI  was  specified  as  the  LM02  TASS  solution 
(SOOkT  surface  burst)  at  /slQmins..  and  the  dispmsion  model  input  was  derived  from  the 
(r4;)-distributions,  as  appropriate.  The  single  burst  release  location  was  SS'N,  37.5*E. 
The  multiburst  scenario  consisted  of  133  idmidcal  clouds,  all  released  at  the  same  initial 
time,  with  the  horizontal  spatial  distribution  shown  in  Figure  2-1.  Three-dimensional 
wind  fields  were  obtained  at  12-hour  intervals  from  the  U.S.  Air  Force  Environmental 
Technical  Applications  Center  (ETAC)  for  16  January  and  4  July  1978.  Tlie  winds  were 
provided  on  a  standard  64  x  64  grid  for  the  Northern  hemisphere,  which  was  then 
interpolated  onto  a  latitude-longitude  mesh  giving  an  effective  ^atial  resolution  of  about 
40Qkm.  Standard  pressure  levels  were  specified  for  the  vertical  profile,  and  a  standard 
U.S.  atmosphere  was  used  to  determine  geometric  heights.  Release  time  was  OOZ  for 
each  of  the  two  days,  and  the  cloud  evolution  was  computed  out  to  16  hours  after  release. 


2.2  SCIPUFF  MODEL  IMPROVEMENTS. 

Two  new  features  were  added  to  SCIPUFF  for  die  'Dust-Off  comparison  exercise. 
An  improved  treatment  of  wind  shear  effects  was  inrarporated,  and  the  capability  to  use  a 
'relative  diffusion'  estimate  was  added  to  the  code.  The  'relative  diffiirion'  prediction 
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Figure  2-1.  Multiburst  scenario  release  locations. 


allowed  a  more  direct  comparison  between  the  SCIPUFF  predictions  and  those  from 
other  models. 

2^1  Inqiroved  Wind  Shear  Treatment 

The  wind  profile  at  the  single  burst  release  location  for  the  winter  case  is  shown  in 
Figure  2-2,  and  clearly  displays  a  very  strong  vertical  shear.  Over  a  16  hour  period,  the 
distortion  due  to  the  vertical  shear  dominates  the  dispersion  of  the  cloud.  The  existing 
version  of  SCIPUFF  used  a  Gaussian  shape  for  each  of  the  collection  of  puffs,  with 
separate  spread  parameters  in  each  of  the  three  orthogonal  directions,  that  is,  East-West, 
Nordi-South,  and  vertical  Vertical  wind  shear  was  represented  as  an  enhanced  horizontal 
spread  rate  for  each  component  direction.  This  method  does  not  properly  describe  wind 
shear,  which  distorts  the  Gaussian  shape  rather  than  diffusing  it,  so  a  more  complete 
treatment  was  implemented. 

The  restriction  of  the  Gaussian  spatial  moment  description  to  the  three  coordinate 
axes  is  equivalent  to  specifying  a  diagonal  form  for  the  second-moment  tensor.  The 
general  definition  of  the  second-moment  is 


4 


Velocity  (m/s) 

Figure  2-2.  Winds  for  16  January  1978  at  55®N,  37.5®E  and  OOZ. 


(2.i: 


where  0  is  the  total  mass  of  the  puff,  c  is  the  local  concentration,  and  Xj  is  the  coordinate 
vector  relative  to  the  puff  centroid.  The  local  puff  concentration  is  represented  by  the 
generalized  Gaussian  distribution 


c(»)« 


&2) 


The  ^neral  Gaussian  distribution  th^fore  requires  6  moments  instead  of  3,  since 
the  second-rank  tensor  is  symmetric.  Tt»  effect  of  wind  shear  on  the  evolution  of  the 
second  moments  is  givra  by 


dt 


(2.3) 


so  that  a  completely  general  wind  field  can  be  included.  The  wind  shear  terms  are  in 
addition  to  the  turbulent  diffusion  terms,  but  it  is  important  to  note  that  a  non-divergmit 
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velocity  field  cannot  increase  the  volume  of  the  puff  in  the  absence  of  diffusion  terms. 
The  volume  is  proportional  to  and  conservation  of  volume  is  an  equivalent 

statement  to  zero-divergence  in  the  velocity  field.  Conservation  is  derivable  from 
Equation  (2.3)  using  the  definition  of  the  determinant 

6Def(<T)  =  (2.4) 

(Jeffieys,  1963).  where  ^  is  the  alternating  third  rank  tensor.  Differentiating  Equation 
(2.4)  with  respect  to  time  gives 

2—{Deti<y))  — 

and  substituting  from  Equation  (2.2),  we  obtain 

d  /  V  {  But  \ 

=  2Dtfr(<T)£,a^£to^  ^ 

^4Det(&)^ 

dXf 

which  is  zero  for  a  divergence-free  velocity  field. 

The  equations  for  the  second-moments  are  advanced  sequentially  by  means  of 
analytical  solutions  for  the  diagonal  strain  components,  followed  by  each  off-diagonal 
velocity  gradient  component  individually.  These  analytic  solutions  conserve  the  puff 
volume  so  the  only  diffusive  effect  is  due  to  the  turbulence  parameterization. 

Generalization  of  the  puff  description  requires  modification  to  the  puff  interaction 
calculation,  local  concentration  computations,  and  also  in  the  puff  splitting/merging 
algorithms.  The  overlap  integrals  and  concratration  evaluations  must  account  for  the 
skewed  axes  of  the  Gaussian,  and  the  splitting  rules  utilize  the  off-diagonal  Gaussian 
moments  to  determine  the  locations  of  the  newly-created  puffs.  Both  splitting  and 
merging  conserves  aU  tte  puff  moments,  and  the  splitting  algorithm  produces  puffs  with 
smaller  volume  and  realizable  moments. 

A  demonstration  of  the  capability  provided  by  the  improved  treatment  of  shear  is 
illustrated  in  Figure  2-3,  which  shows  the  single  burst  cloud  evolving  over  8  hours  in  the 
January  wind  field.  The  perspective  views  of  the  10~^*  gcm~^  iso-surface  use  a  distorted 
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scale  to  represent  the  horizontal  stretching  of  the  cloud,  which  covers  a  range  of  several 
hundred  kUometers  by  the  end  of  the  computation.  The  ability  of  the  puffs  to  distort 
under  the  action  of  the  persistent  shear  allows  the  code  to  maintain  a  continuous 
representation  of  the  cloud  without  unreasonable  computational  effort  This  calculation 
contained  646  puffs  at  the  last  time  illustrated  with  10  paiticle  size  groups  in  the 
representation. 

122  Rdative  IMffiision  PredictioiL 

The  standard  concentration  field  output  from  SCIPUFF  consists  of  an  ensemble 
average  concentration  along  with  a  prediction  of  the  variance  about  the  mean.  The 
ensemble  is  depmident  on  the  input  wind  field  data,  which  is  always  incomplete  due  to 
limited  spatial  and  temporal  resolution,  and  the  model  predictions  give  a  probabilistic 
description  of  the  concentration  field  resulting  from  the  uncertain  winds.  This  makes 
comparison  between  SCIPUFF  predictions  and  those  from  a  conventional  deterministic 
model  very  difficult  In  order  to  facilitate  such  comparisons,  we  modified  the  SCIPUFF 
diffusion  algorithm  to  allow  the  use  of  a  'relative  diffusion'  scheme.  The  concepts  of 
'absolute'  and  'relative'  dispersion  are  discussed  by  Pasquill  (1974)  and  can  loosely  be 
thought  of  as  continuous  'plume'  dispersion  and  instantaneous  'puff  dispersion, 
respectively.  Absolute  dispersion  corresponds  to  the  normal  SCIPUFF  ensemble 
average,  where  (unresolved)  eddies  of  all  scales  are  included  in  the  dispersion  process. 
This  can  be  interpreted  as  a  long  time  average  of  a  continuous  plume,  where  the 
averaging  is  sufficiently  long  to  sample  all  the  eddies.  More  generally,  the  ensemble 
average  is  a  probabilistic  prediction,  where  large  eddies  are  regarded  as  meandering  the 
cloud  and  introducing  uncertainty  into  its  position;  a  long  term  average  will  reduce  this 
uncmtainty  and  give  the  usual  plume  average,  but  a  short  term  average  can  still  contains 
significant  uncertainty. 

We  can  regard  a  deterministic  prediction  of  the  cloud  dispersion  as  an  attempt  to 
predict  the  likely  concentration  in  the  cloud,  ignoring  the  uncertainty  in  its  position. 
There  is  no  clear  distinction  between  the  'meandering'  eddies  and  the  'diffusing'  eddies, 
but  there  is  a  natural  way  to  restrict  the  SCIPUFF  prediction  to  include  only  the  cloud- 
scale  eddies.  SCIPUFF  already  requires  an  estimate  of  the  cloud-scale  turbulence  to 
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a)  f  =  0  hrs 


Hguie  2-3.  SCIPUFF  results  for  the  single  burst  scenario  using  the  January  wind  field. 
Perspective  views  of  the  lO"'®  g/cm^  iso-surface  at  r  =  0,1  and  8  hrs. 
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calculate  the  dissipation  of  die  concentratioo  fluctuation  variance,  so  we  define  a  'relative' 
diffusion  prediction  to  utilize  only  this  part  of  the  turbulence  spectrum.  Specifically, 
given  a  cloud  length  scale,  Xc ,  and  a  turbulence  length  scale  representing  the  'maximum' 
energy-containing  eddy  size.  Ah  ,  we  represent  the  velocity  variance  available  for 
diffusion  as 


(2.6) 


for  Ac  <  Ah  .  whereug  is  the  ambient  velocity  variance.  This  representation  is  consistent 
with  the  observed  S/3-power  law  spectrum  (Nastrom  and  Gage,  1985). 


23  TEST  CASE  RESULTS. 

Results  from  the  model  inter-comparison  exercise  are  not  the  subject  of  this 
report,  but  we  shall  present  several  examples  of  the  SCIPUFF  output  for  the  designated 
test  cases.  As  discussed  above,  single  and  multiple  burst  scenarios  were  specified  for  two 
meteorological  conditions.  Some  experimentation  was  performed  to  determine  the 
appropriate  resolution  for  the  calculation,  particularly  with  regard  to  the  particle  size 
distribution.  A  logarithmically  uniform  size  bin  spacing  was  employed,  ranging  from 
1pm  up  to  l()24pm,  with  a  dust  density  of  2.0gcm~^.  Test  runs  were  made  with  5, 10, 20, 
and  40  size  bins,  showing  that  10  bins  was  sufficient  to  maintain  an  accurate  description 
of  the  cloud  over  the  time  period  of  the  calculation.  The  model  resolution  was  effectively 
extmided  by  the  inclusion  of  an  additional  vertical  spread  term  to  represent  the  range  of 
particle  fall  velocities  over  the  finite  range  of  the  aze  bin.  Thus,  the  particles  fall  with 
the  mass-weighted  mean  velocity,  given  by 
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vj  =  2^a,v^- 

ial 

where  Vdi  =  v^dmin) 

Vd2  —  Vd(dmeaa) 

V<o  =  Vd(daiMJi) 

referring  to  the  particle  fall  speed  at  the  minimum,  mean,  and  maximum  sizes  in  the  size 
bin,  respectively.  The  mass  weights,  Oi ,  are  inversely  proportional  to  the  diameter, 
consistent  with  the  -4  power  law  for  the  number  density.  The  vertical  puff  spread 
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equation  includes  a  growth  term,  so  that  Oi  increases  like  ,  where  the  r.in.s.  fall 
velocity  is  obtained  from  the  three  values  at  the  bin  boundaries  and  bin  center.  Le, 

i=l 

This  additional  tmm  spreads  the  puffs  and  maintains  a  continuous  vertical  concentration 
distribution,  reducing  the  tendency  for  puffs  to  separate  due  to  differential  fall  speeds. 

Horizontal  cross-sections  through  the  concentration  field  at  tsShrs  and  at  three 
standard  altitudes  ate  presented  in  Figures  2-4  through  2-7.  These  fields  show  both  total 
dust  concentration  and  also  the  concentration  of  dust  particles  with  sizes  greater  than 
SOpm  from  each  of  the  four  cases,  Le.,  winter  and  sununer  winds,  single  and  multiburst 
scenarios.  These  four  calculations  were  made  using  the  relative  diffusion  algorithm 
described  in  Section  2.2.2,  and  there  are  several  features  of  interest  in  the  results. 

Hrst,  we  note  the  continuous  nature  of  the  SQPUFF  prediction.  The  particle  size 
distribution  resolution  together  with  the  generalized  Gaussian  puff  description  allows  the 
code  to  maintain  an  accurate  description  of  the  highly  distorted  cloud  as  it  shears  and 
separates  the  different  particle  sizes  This  is  particularly  evident  in  the  single  burst, 
winter  wind  case  shown  in  Hgure  2-4.  Here  the  strong  wind  shear  stretches  the  cloud 
over  a  distance  of  order  1000km  over  the  8hr  period.  This  stretching  sqtarates  the  main 
cloud  fiom  the  pedestal,  as  can  be  seen  by  comparing  the  total  dust  plot  with  that  for  the 
dost  above  50p.m.  The  larger  particles  in  the  section  at  z=100m  must  necessarily  have 
originated  at  much  higher  altitudes,  since  they  have  been  falling  for  8hrs,  while  the 
smaller  particles  originated  at  altitudes  much  closer  to  lOQm.  Comparison  between  the 
two  concentration  fields  shows  that  the  larger  particles  form  the  Western  section  of  the 
cloud  furthest  from  the  burst  location,  while  the  smaller  pedestal  particles  make  up  the 
high  concentrations  in  the  North-East  These  cross-section  plots  can  be  compared  with 
the  perspective  view  of  the  cloud  in  Figure  2-3,  where  the  "curtain"  of  larger  particles 
falling  out  of  the  sheared  cloud  downstream  of  the  pedestal  is  clear. 

The  summer  wind  case  in  Hgure  2-5  shows  a  simpler  cloud  structure,  but  the 
wind  shear  effect  is  still  obvious.  The  multiburst  scenarios  (Figures  2-6  and  2-7)  can  be 
thought  of  as  a  simple  superposition  of  133  single  burst  clouds  with  appropriate  origin 
shifts.  These  calculations  ^ow  tiie  same  basic  features  as  the  single  burst  scenarios,  but 
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a)  Dust  larger  than  50  pm 


b)  Total  dust 
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Hgttie  2-4.  SCIPUFF  results  for  the  ^gle  burst  scenario  using  the  January  wind  field. 
Horizcmtal  cross-sections  at  r=  8  hrs.  Contour  levels  of  10'^^  lO*^^,  lO'^L 
10-*®,  1(J®,  10^  and  10-^  gm/cm^. 


there  is  less  detail  in  the  clouds  since  the  oveiiapping  concentration  fields  trad  to  smooth 
out  the  individual  features. 

The  intercepted  dust  along  an  aircraft  flightpath  is  an  important  quantity  for 
strategic  planning  purposes,  and  Figure  2-8  shows  the  SCIPUFF  prediction  for  both 
Nmth-South  and  East-West  paths  through  the  single  burst  cloud  shown  in  Figure  2-4. 
Integrated  dust  concentrations  are  shown  at  the  three  altitudes  as  a  fraction  of  intercept 
latitude  or  longitude.  The  shear-induced  offset  is  clear  in  the  altitude  dependence,  and  the 
near-surfBoe  separation  between  tire  main  cloud  and  the  pedestal  is  also  evictent,  with  the 
pedestal  dust  lying  to  the  North  and  the  East  of  the  main  cloud. 
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a)  Oust  larger  than  50 


b)  To^l  dust 


Longitude 


Hgure  2-5.  SCIPUFF  results  for  the  single  burst  scenario  using  the  July  wind  field. 

Horizontal  cross-sections  at  r  =  8  hrs.  Contour  levels  of  lO'^^  lO-^^,  10-“, 
10-^0, 10^,  10^  and  lO-^gm/cm^. 


As  noted  earlier,  the  SCIPUFF  predictions  shown  above  were  obtained  using  the 
relative  diffusion  model  to  facilitate  comparison  with  deterministic  predictions. 
However,  the  normal  SCIPUFF  prediction  is  an  absolute  diffusion  estimate  of  the  mean 
concentration  together  with  a  calculation  of  the  expected  variance  about  this  mean.  The 
absolute  diffusion  prediction  shows  a  much  wido*  spread,  since  all  scales  of  eddies  ate 
included  in  the  dispersion  process  which  now  represents  positional  uncertainty  in  addition 
to  cloud  diffusion.  Figure  2-9  shows  the  z=100m  fields  from  the  single  burst  test  case 
with  winter  winds.  The  mean  concentration  contours  in  this  figure  can  be  contrasted  with 
those  in  Figure  2-4.  The  increased  diffusion  is  obvious,  but  Figure  2-9  also  contains  the 
rjn.s.  fluctuation  prediction  and  this  is  seen  to  be  larger  than  the  mean  value,  in  general 
This  indicates  a  large  uncertainty  in  the  predicted  concentration  at  any  particluar  location 
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a)  Dust  larger  than  50  tun 
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Figure  2-6.  SQPUFF  results  for  the  multi-burst  scenario  using  the  January  wind  field. 

Horizontal  cross-sections  at  r  =  8  hrs.  Contour  levels  of  10-^^  10-‘^ 

10-1®,  10-®,  10^  and  lO-”^  gm/cm^. 


45 


and  should  be  interpreted  as  mostly  due  to  uncertainty  in  position,  i.e.,  the  available  wind 
information  is  not  sufficient  to  reliably  predict  the  cloud  position  relative  to  the  actual 
cloud  widdi. 

The  fluctuation  intensity  characterizes  the  random  nature  of  the  cloud,  and  this  can  be 
used  to  generate  'realizations'  of  the  instantaneous  concentration  distribution.  We  use  the 
statistical  description  together  with  a  'fractal'  representation  technique  to  generate  a 
random  field  with  an  appropriate  spatial  structure.  The  random  fields  created  by  this 
technique  possess  the  same  statistical  properties  as  the  SCIPUFF  prediction  but  other 
details  must  be  assumed.  Thus,  we  use  a  clipped  normal  representation  of  the  probability 


13 


a)  Dust  larger  than  50  pm 
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Figure  2-7.  SCIPUFF  results  for  the  multi-burst  scenario  using  the  July  wind  Held. 

Horizontal  cross-sections  at  r  =  8  hrs.  Contour  levels  of  10-^^  10-^^ 

10-‘0  10-9  icHand  IQ-^gm/cm^. 


distribution  function  for  concentration  at  any  point,  together  with  a  spatial  correlation 
scale  based  on  the  SCIPUFF  prediction  of  the  cloud  scale,  A^,  and  we  also  specify  the 
fractal  dimension  of  the  field  so  that  a  concentration  contour  has  dimension  1.3.  The 
fractal  dimension  is  based  on  estimates  firom  atmospheric  observations  of  natural  clouds 
(Lovejoy,  1982).  The  fractal  fields  are  intended  to  represent  possible  clouds  from  the 
ensCTible  of  aU  possible  clouds,  and  two  examples  are  given  in  Figure  2-10.  These  fields 
were  geirerated  from  the  statistics  shown  in  Hgure  2-9  but  are  clearly  very  different  firom 
the  average  concentration  field.  The  intermittent  nature  of  the  clipped  normal  distribution 
gives  much  smaller  regions  of  high  concentration.  The  concentration  levels  and  areal 
coverage  are  broadly  similar  to  the  relative  diffusion  prediction  in  Figure  2-4,  but 
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Figure  2-8.  Intercepted  dust  masses  along  horizontal  flight  paths  through  the 
single  burst  cloud  after  8  hours  with  winter  winds,  (a)  East-West 


paths,  (b)  North-South  paths. -  zslOOm; 
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Figure  2-9.  SCIPUFF  results  for  the  single  burst  scenario  using  the  January  wind  field 
and  absolute  diffusion  scheme.  Horizontal  cross-sections  atr  =  8  hrs  and 
z  =100  m.  (a)  mean  total  dust  concentration,  0>)  r.m.s.  fluctuation  in 
total  dust.  Contour  levels  of  10-^  10*^2^  lO-^L  10*^®  and  lO*®  gm/cm^ 


Figure  2-10.  Fractal  realizations  from  the  ensemble  statistics  of  Figure  2-9. 

Horizontal  cross-sections  at  rs  8  hrs  andz^lOO  m.  Contour  levels 
of  10-13  10-12  10-11, 10-iOand  10^  gm/cml 
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the  Helds  are  much  more  complex  due  to  the  random  nature  of  the  variability. 
Unfortunately,  the  intermittency  of  the  clipped  normal  distribution  gives  very  sharp 
boundaries  around  the  non-zero  concentration  regions,  making  it  difficult  to  distinguish 
the  logarithmic  contour  levels.  However,  the  maximum  concentration  values  in  both 
realizations  are  close  to  KHg/cc  and  occur  in  the  region  around  52*N  36*E,  with  levels 
roughly  an  order  of  magnitude  smaller  in  the  South-Western  regions  of  the  cloud. 
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SECTION  3 


'ANATEX'  MODEL  EVALUATION  STUDY 

Model  validation  is  a  vital  part  of  any  program  development  but  this  has  always 
been  a  difficult  exercise  for  long  range  atmospheric  dispersion  models.  Large  scale 
atmospheric  flow  fields,  with  length  scales  ranging  from  hundreds  to  thousands  of 
kilometers,  cannot  be  simulated  in  the  laboratory  with  sufficient  verisimilitude  for 
dispersion  effects  to  be  measured.  We  must  therefore  rely  on  field  observations  for 
model  evaluation  purposes.  Field  programs  must  deal  with  the  inherent  uncertainties  of 
the  wind  field,  which  carmot  be  measured  with  enough  detail  for  reliable  trajectory 
calculations,  and  also  with  limited  sampling  of  the  tracer  distribution  over  the  large 
horizontal  scales  of  interest  Most  model  evaluations  are  therefore  restricted  to 
qualitative  comparison  at  continental  scales. 

The  Across  North  America  Tracer  Experiment  (ANATEX)  was  chosen  to 
evaluate  the  SCIPUFF  model  predictions  because  it  is  an  unusually  complete  study  of 
long-range  transport  and  di^)^oa  The  wind  field  uncertainties  are  unavoidable,  but  the 
ANATEX  data  provides  an  exceptionally  large  data  base  for  model  evaluation  at  the 
scales  of  interest  here.  This  is  of  particular  valt^  in  evaluating  a  probabilistic  di^rsion 
model  like  SCIPUFF,  where  an  assessment  of  the  uncertainty  in  the  prediction  is  made, 
since  the  comparison  with  observations  can  only  be  made  on  a  statistical  basis. 
ANATEX  was  conducted  between  January  and  March,  1987,  under  the  auspices  of 
NOAA  (National  Ocean  and  Atmosphere  Administration).  Specialized  tracer  materials 
were  released  from  two  sources  in  Montana  and  Miimesota  in  cycles  of  short  duration 
over  the  pmod  of  the  experiment  The  tracers  were  detected  on  a  networit  of  77  surface 
samplers,  illustrated  in  Figure  3-1,  covering  the  entire  Eastern  half  of  the  United  States 
and  the  Canadian  Maritime  Provinces.  The  samplers  recorded  24-hour  average 
concentrations  throughout  the  experiment  In  addition,  three  aircraft  made  direct 
sampling  flights  through  the  plumes  within  a  few  hundred  kilometers  of  the  source  on  a 
number  of  occasions. 

Wind  fields  were  provided  by  NOAA  firom  the  Nested  (3tid  Model  (NGM),  an 
operational  forecast  model  with  a  spatial  resolution  of  about  l(X)km  over  the  continental 
United  States.  Three-dimensional  wind  fields  were  stored  every  2  hours  throughout  the 
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Hgure  3-1.  ANATEX  release  and  surface  sampler  locations. 


experimental  period,  with  the  exception  of  one  week  in  March  where  the  meteorological 
data  is  incomplete.  The  data  was  therefore  divided  into  two  parts,  covering  7  weeks  of 
January  and  February,  and  3  weeks  of  March,  where  wind  fields  were  available. 
SCEPUFF  was  run  continuously  for  these  periods,  simulating  the  multiple  releases  from 
each  of  the  two  release  locations,  and  tracking  the  total  tracer  concentrations  over  the 
Eastern  U.S.  A  complete  time  history  of  the  tracer  concentrations  at  each  surface  sampler 
site  was  saved,  and  24-hour  averages  were  computed  for  comparison  with  the  surface 
observations.  In  addition,  tracer  concentrations  along  the  specified  aircraft  fligh^ths 
were  also  saved  for  comparison  with  the  direct  sampling  by  the  aircraft  systems.  The 
complete  results  from  the  validation  study  were  published  in  the  Journal  of  Applied 
Meteorology,  and  the  paper  is  included  as  an  Appendix  to  this  report.  The  results  will 
therefore  only  be  summarized  in  this  section,  and  the  reader  is  referred  to  the  Appendix 
for  further  details. 
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A  Q'pical  comparison  between  the  fnedicted  surface  concentration  distribution  and 
the  observed  values  is  shown  in  Figure  3-2.  The  general  shape  and  magnitude  of  the 
surface  concentrations  is  correctly  predicted  by  the  model,  but  detailed  comparison  of  the 
concentrations  at  specific  sampler  locations  is  not  very  good.  The  reason  for  the  poor 
correlation  is  mostly  the  inaccuracy  of  the  cloud  trajectory  prediction.  Small  errors  in  the 
cloud  position  produce  large  discrepancies  between  the  predicted  and  observed 
concentrations,  even  though  the  overall  distribution  may  be  very  similar.  This  is  a 
common  problem  in  model  evaluation,  and  various  techniques  have  been  developed  to 
remove  the  effect  of  trajectory  errors.  One  of  the  standard  comparison  statistics  is  the 
cumulative  concentration  distribution,  i.e.,  the  function  F{c)  defined  such  that  F  is  the 
fraction  of  the  observations  less  than  the  value,  c  .  The  function  can  be  defined  for  any 
subset  of  the  observations  to  examine  dependence  on  location  or  time,  but  the  comparison 
becomes  more  meaningful  for  larger  samples  since  the  statistical  uncertainties  are 
reduced. 

A  comparison  between  the  SCIPUFF  predictions  of  the  24-hour  average  tracer 
concentrations  over  the  whole  network  for  the  March  period  is  shown  in  Figure  3-3.  The 
agreement  is  seen  to  be  good  over  the  entire  range  of  concentrations,  indicating  that  the 
correct  fraction  of  concentrations  is  predicted  at  each  level.  These  predictions  were 
obtained  using  the  boundary  layer  description  provided  by  NOAA,  but  results  using  this 
estimate  were  not  as  good  for  the  January-February  period.  The  NOAA  estimates  of  the 
boundary  layer  depth  appeared  to  be  very  low  for  this  period,  and  a  SCIPUFF  prediction 
using  a  standard  1000m  deep  layer  during  the  daytime  gave  much  improved  agreement 
The  comparison  for  individual  sampler  arcs  is  shown  if  Figure  3-4,  illustrating  that  the 
range  dependence  of  the  concentration  values  is  also  correctly  predicted.  Further 
discussion  of  these  results  can  be  found  in  the  Appendix. 

The  uncertainty  in  the  24-hour  samples  is  mainly  due  to  the  trajectory  differences 
arising  from  the  numerical  forecast  errors.  We  have  not  attempted  to  model  the  random 
forecast  errors,  but  have  simply  compared  the  cumulative  distributions.  The  contribution 
from  unresolved  eddies  (smaller  than  the  ICXlkm  NGM  grid)  is  insignificant  in  the  surface 
tracer  data  since  the  24-hour  average  is  long  enough  to  reduce  these  uncertainties  to 
negligible  proportions.  This  is  not  true  for  the  short-term  aircraft  sampling,  which  is  only 
a  6-minute  duration.  The  uncertainty  in  these  values  was  therefore  a  significant  factor  in 
the  model  evaluation.  A  direct  comparison  between  the  predicted  and  observed 
distributions  showed  a  large  overestimate  of  the  number  of  trace  (less  than  lfl/1) 
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Hgure  3>2.  Comparison  of  the  SCIPUFF  24hr  average  PTCH  surface  concentration  patterns 
with  the  observed  patterns  for  the  period  23-25  March. 
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Hgure  3-3.  Predicted  and  observed  ground  sample  cumulative  distribution  functions  for 
die  period  9-26  March. 
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Hguie  3-4.  Predicted  and  observed  ground  sample  cumulative  distribution  functions  for 
PTCH  along  fliree  arcs  for  the  p^od  5  January-21  February.  Predictions 
were  made  using  a  standard  boundary  layer  depth  of  lOOQm. 
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concentratioii  samples,  as  can  be  seen  in  Figure  3-5.  The  predicted  confidence  bounds 
are  included  in  the  Hgure,  but  it  is  clear  that  the  observations  lie  well  outside  the 
acceptable  range.  A  likely  reason  for  the  discrepancy  is  the  non-random  nature  of  the 
aircraft  sampling;  flights  were  made  on  the  basis  of  forecast  trajectories,  with 
readjustment  after  initial  samples  were  taken.  We  expea  the  aircraft  to  record  fewer  zero 
concentrations  than  a  completely  random  sample,  therefore,  since  there  was  a  definite 
effort  to  locate  the  actual  plume.  We  can  only  compare  the  distribution  of  non-zero 
concentrations,  in  this  case,  since  the  conditional  probability  of  intercepting  the  plume 
caimot  be  accurately  estimated.  If  we  adjust  the  predicted  zero  fraction  to  match  the 
observations  we  obtain  the  improved  agreement  shown  in  Figure  3-6,  which  is  much 
more  acceptable.  Most  of  the  observed  distribution  now  falls  within  the  predicted 
confidence  bounds  for  the  sample. 
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Figure  3-6.  Adjusted  predicted  and  obsoved  aircraft  sample  cumulative  distribution 
functions  for  the  period  9-28  March. 

As  part  of  the  ANATEX  study,  we  also  examined  the  possibility  of  predicting  the 
expected  distribution  of  concentrations  based  simply  on  the  long  term  average  wind 
statistics  for  the  entire  Eastern  United  States.  There  is  no  detailed  information  on  the 
wind  field,  only  average  speed  and  direction  as  a  function  of  altitude,  together  with  a 
variance  about  the  mean.  The  SCIPUFF  prediction  is  compared  with  the  surface  sampler 
observations  in  Hgure  3-7,  and  shows  reasonable  agreemoit  considering  the  very  limited 
information  used  as  input  The  model  predicts  a  large  uncertainty  range,  due  to  the  very 
uncertain  wind  field  input  so  the  expected  concentration  distribution  is  sampled  from  the 
predicted  range  with  an  assumed  probability  density  function.  For  all  our  shorter  range 
(less  than  l(X)km)  predictions  with  shorter  averaging  times  (up  to  1  hour),  we  have  found 
that  the  'clipped-normal'  distribution  gives  a  good  prediction,  but  the  results  in  Hgure  3-7 
were  obtained  using  a  log-normal  distribution.  The  clipped-normal  gave  too  much 
intermittency  and  was  unable  to  match  the  observations.  There  is  some  indication, 
therefore,  that  the  longer  averaging  and  transport  times  are  mors  appropriately 
represented  by  a  log-normal  probability  distribution.  The  results  do  indicate,  however, 
that  reasonable  estimates  of  expected  concentration  values  can  be  predicted  with  limited 
descriptions  of  the  wind  field.  Particular  events  cannot  be  predicted  without  detailed 
information,  but  the  correct  probability  of  a  particular  event  can  be  predicted. 
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Figme  3-7.  Predicted  and  obs^ed  ground  sample  cumulative  distribution  functions  for 
tbe  period  9-28  March.  Predictions  were  made  using  the  ensemble 
meteorology. 
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SECTION  4 


BOUNDARY  LAYER  DIFFUSION  EFFECTS 

4.1  BACKGROUND. 

The  blast  of  a  large  nuclear  weapon  lofts  large  amounts  of  dust  and  ejecta  high 
into  the  atmo^here.  In  addition  to  this  high  cloud  of  material,  the  strong  winds  near  the 
surface,  induced  by  the  rapidly-rising  nreball,  scour  loose  material  to  form  a  broad 
pedestal  of  low-lying  dust  The  flow  associated  with  the  blast  of  the  explosion  also  lofts 
the  sand/soil  particles,  so  that  the  pedestal  can  represent  a  significant  fraction  of  the  total 
mass  of  material  injected  into  the  atmosphere  by  the  burst  The  pedestal  material  is 
confined  to  a  relatively  fallow  layn*  above  the  ground,  but  the  action  of  turbulent  eddies 
on  the  dust  can  lift  it  to  heights  of  order  1km,  which  represents  a  typical  mixing  depth  for 
surface-driven  turbulence.  The  possibility  that  the  large  dust  mass  in  the  pedestal  can 
remain  within  the  surface  mixing  layer  for  periods  of  hours  raises  concern,  and  a  mote 
detailed  study  than  has  hitherto  been  undertaken  is  needed  to  answer  questions  about  the 
importance  of  this  phenomenon. 

There  are  two  basic  physical  mechanisms  determining  the  evolution  of  the  dust 
field  in  the  lowest  layers  of  the  atmosphere;  first,  the  turbulent  eddies  generated  by  heat 
input  from  the  ground  and  by  mechanical  drag  on  the  rough  surface  transport  the  particles 
through  the  air,  and  second,  the  deposition  of  the  particles  onto  the  ground.  The  first 
mechanism  tends  to  maintain  the  particles  in  suspension,  while  the  second  is  purely  a 
reduction  of  airborne  material.  The  relative  balance  between  these  two  mechanisms 
depends  strongly  on  the  particle  size.  For  large  particles,  i.e.  larger  than  100p.m, 
gravitational  settling  is  the  major  deposition  mechanism,  and  the  turbulent  eddy  velodty 
must  be  compared  with  the  equilibrium  fall  rate  to  determine  the  importance  of  turbulent 
convection  on  the  paiUcles.  Typical  eddy  velocities  are  of  the  order  of  O.Sms~^  while 
particles  of  100pm  size  have  a  terminal  velocity  of  about  0.3ms~^  so  that  particles  larger 
than  100pm  are  unlikely  to  be  prevented  from  falling  for  very  long.  This  preliminary 
balance  determines  the  size  range  of  particles  which  can  be  influenced  significantly  by 
atmospheric  boundary  layer  turbulence  to  be  below  100pm.  For  smaller  particles  the 
gravitational  settling  rate  decreases  rapidly  with  diameter,  but  under  turbulent  conditions 
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the  most  important  deposition  mechanism  can  be  tte  impact  of  particles  onto  the  surface 
roughness  elements.  Thus,  a  realistic  examination  of  die  response  of  the  dust  cloud  under 
turbulent  conditions  needs  to  account  for  the  turbulent  deposition  of  the  smaller  particles. 

In  this  study,  we  focus  on  the  vertical  distribution  of  the  dust  mass,  i.e.  we  are  not 
concerned  with  the  horizontal  transport  and  dispersion  of  the  cloud.  The  conservation 
equation  for  the  dust  concentration  can  be  integrated  over  the  entire  horizontal  plane  to 
eliminate  the  horizontal  terms,  and  reduce  the  problem  to  a  one-dimensional,  time- 
dependent  set  of  equations.  The  wind  transport  and  horizontal  di^iision  only  redistribute 
the  mass  in  the  horizontal,  so  that  we  can  ignore  these  effects  if  we  assume  that  the 
boundary  layer  is  horizontally  homogeneous  and  we  are  only  interested  in  the  vertical 
distribution.  The  one-dimensional  results  from  this  study  will  therefore  give  information 
on  how  fast  the  dust  is  removed  from  the  atmosphere  in  an  overall  sense,  and  on  the 
altitude  dependence,  but  provides  no  guidance  in  assessing  an  actual  local  concentration 
level  since  this  will  depend  on  the  horizontal  redistribution. 

The  next  section  briefly  describes  the  model  equations  used  to  represent  the 
atmospheric  turbulence  and  the  surface  deposition  rate,  together  with  the  numerical 
techniques  employed  in  their  solution.  Section  4.3  presents  the  results  of  the  calculations, 
which  cover  a  range  of  conditions,  including  wind  speed,  surface  heating,  surface 
roughness,  time  of  release,  and  particle  size  spectrum.  Finally,  Section  4.4  discusses  the 
results,  and  presents  the  general  conclusions  from  the  investigation. 


4,2  NUMERICAL  MODEL. 

4.2.1  Turbulent  Transport  Model. 

The  model  used  in  this  study  is  based  on  turbulent  dispersion  research  in  the  field 
of  particulate  transport  from  industrial  sources.  The  basic  turbulence  model  is  derived 
from  the  second-order  closure  technique  and  a  surface  deposition  formulation  consistent 
with  the  turbulence  closure  is  employed.  Both  models  are  fully  described  by  Lewellen 
and  Sheng  (1980). 

The  particle  size  spectrum  is  represented  by  a  number  of  discrete  bins  with  no 
interaction  between  different  size  groups,  so  that  the  bins  can  be  treated  independently. 
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This  neglects  any  particle  agglomeration  or  break-up,  although  tlmse  mechanisms  could 
conceivably  be  significant  in  the  high  density  environment  of  the  early  cloud.  However, 
a  description  of  these  effects  is  beyond  the  current  state-of-the-art  and  is  therefore 
omitted  from  the  present  study.  Similar  remarks  apply  to  electrostatic  effects,  which 
could  be  important  in  the  transport  of  very  small  particles,  but  very  little  is  known  about 
their  role  in  the  atmosphere. 

We  denote  the  (horizontally-integrated)  concentration  in  a  size  group 
characterized  by  diameter,  da,  by  c^a)  >  tbe  conservation  equation  for  C(a)  as 

dt  dz  dz 

where  the  overbar  denotes  the  Reynolds  ( or  ensemble)  average  operator,  and  the  prime 
denotes  the  random  turbulent  fluctuation  about  the  average.  The  second-order  correlation 
term  is  the  turbulent  flux  of  concentration,  which  is  often  parameterized  as  a  diffusive 
flux  proportional  to  the  mean  gradient  in  simple  models.  This  assumption  is  iiuiccurate  in 
some  circumstances,  especially  in  buoyancy-driven  flows,  and  at  best  requires 
specification  of  an  eddy  diffusivity.  The  second-order  closure  technique  derives  a 
predictive  equation  for  the  turbulent  flux  itself,  and  seeks  to  rationaUy  model  the  higher- 
order  moments  which  iqrpear  in  this  equation. 

The  gravitational  settling  velocity  V(a)  ,  appropriate  for  the  particle  size,  is 
determined  from  the  balance  between  aerodynamic  drag  and  gravitational  forces,  i.e. 

whm  Pd  is  the  dust  density,  is  the  air  density,  V  is  the  particle  volume,  Ap  is  the 
particle  cross-section  area,  and  cd  is  the  drag  coefficient  The  drag  coefficient  is  a 
function  of  Reynolds  number,  ^e^dlv  ,where  d  is  the  effective  particle  diameter  and  v 
is  the  kinematic  viscosity  of  air,  given  by 

c^=ii(l+0.15Re®")  (4.3) 

Re 

for  Re<l()00,  and  held  at  a  constant  value  of 0.424  for  all  larger  Reynolds  numbers. 

The  equation  for  the  vertical  concentration  flux  is  written  in  the  form 
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(4.4) 


where  9  is  the  rms  turbulent  velocity  fluctuation,  A  is  the  turbulence  length  scale,  g  is  the 
gravitational  acceleration,  Tq  is  the  reference  temperature  for  the  Boussinesq 
approximation,  and  0  '  is  the  temperature  fluctuation.  A  and  Vc  are  dimensionless 
empirical  constants  used  in  the  turt>ulence  modeling,  and  take  the  values  0.75  and  0.3 
respectively.  It  can  be  seen  that  the  boundary  layer  turbulent  velocity  fluctuations  ^pear 
directly  in  tl»  concentration  flux  equatioiL  The  boundary  layer  dynamics  are  computed 
along  with  the  dust  equations,  so  that  the  diurnal  variation  of  the  turbulence  is  correctly 
included  in  the  calculatiotL  The  details  of  the  diurnal  meteorology  as  specified  for  this 
study  will  be  givmi  below. 


AJ12  Turbulent  Deposition  Model. 


The  deposition  of  small  particles  onto  the  surface  can  be  modeled  using  a 
turbulent  deposition  velocity,  Vd  ,  which  augments  the  gravitational  deposition  velocity, 
V(a).  The  additional  flux  of  concentration  onto  the  surface  due  to  turbulent  transfer  is  then 
written  as  the  product  of  the  near-surface  concentration  and  the  turbulent  deposition 
velocity,  which  depends  on  the  character  of  the  surface,  the  particle  size,  and  also  on  the 
turbulence.  We  use  a  relatively  simple  description  of  the  surface  in  terms  of  a  single 
roughness  length  in  this  study;  e.g.  a  forest  is  represented  by  a  large  toughness  length 
compared  with  a  smooth  surface  such  as  desert  sand.  The  large  roughness  length  not 
only  increases  the  surface  drag,  and  hence  the  turbulence  intensity,  but  also  produces  a 
m(»e  efficient  scavenging  of  the  small  particles  as  they  are  swept  through  the  roughness 
elemmits  ( i.e.  leaves  in  the  forest  case). 


The  specific  form  of  the  turbulent  deposition  velocity  is  obtained  from  the  work 
of  Lewellen  and  Sheng  (1980).  We  base  the  formulation  on  the  passive  scalar  flux  rate, 
modified  as  tq)propriate  for  finite  particle  size.  Thus  we  write 


(4.5) 


where  u  is  the  mean  wind  at  a  reference  height  ( the  lowest  model  grid-point  in  our 
numerical  integration),  ce  is  the  bulk  transfer  coefficient  for  a  passive  scalar,  and  Pj 
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iqxtsents  the  particle  size  effects,  ca  depends  on  the  surface  roughness,  and  the  Monin* 
Obhukhov  length,  which  is  a  measure  of  the  buoyancy  forces  in  the  turbulent  transfer 
process. 


The  particle  effects  are  modeled  as 


P,=0.8 


{DIvf 


(4.6) 


where  Z>  is  the  particle  diffiisivity  due  to  Brownian  motion  (from  Ftiedland^,  1977),  v  is 
the  molecular  viscosity  of  air,  9  is  the  rms  turbulent  velocity,  u*  is  the  surface  friction 
velocity,  and  A  is  the  turbulence  length  scale,  ty  is  the  particle  response  timescale  and  is 
obtained  from  the  drag  law,  so  that 


T  -  JPifL 


(4.7) 


udiere  the  Reynolds  number  is  based  on  q.  The  turbulent  velocity  scale  is  used  in  this 
definition  ance  we  are  primarily  interested  in  modeling  the  response  of  the  particle  to  the 
turbulent  eddy  field.  We  note  that  for  small  particles  with  small  Reynolds  numbers,  the 
timescale  is  independent  of  q  and  depends  only  on  the  fluid  properties  and  the  particle 
size  and  density. 


The  physics  of  the  turbulent  particle  deposition  model  are  governed  by  the 
balance  between  the  inertial  effects  of  the  particles  in  the  accelerating  turbulent  field  and 
Brownian  motion  effects.  Very  small  particles  are  diffused  across  the  viscous  sublayer 
onto  the  surface,  while  larger  particles  are  impacted  initially  onto  the  surface  elements. 
Examples  of  the  variation  of  turbulent  deposition  velocity  with  particle  size  are  shown  in 
Figure  4-1,  and  compared  with  the  gravitational  settling  velocity.  These  results  are 
computed  for  a  reference  height  of  IQm,  and  surface  friction  velocity  of  0.4ms~i  and  a 
particle  density  of  2g/cc.  Two  surface  roughness  values  are  shown,  and  it  can  be  seen 
that  the  turbulent  deposition  only  slightly  exceeds  the  gravitational  term  for  the  low 
roughness  case,  and  this  occurs  only  for  particle  diameters  between  about  Spm  and 
2()pm.  For  the  high  roughness  case,  however,  the  turbulent  deposition  is  almost  a  factor 
of  ten  larger  than  gravitational  settling  for  lOpm  particles,  and  is  the  principal  deposition 
mechanism  for  sizes  between  24m  and  40iim.  The  deposition  velocities  are  of  the  order 
of  a  few  cms~^  for  particles  above  10pm  in  \b&  high  roughness  case.  These  vdocities 
increase  witii  increasing  friction  velocity,  i.e  with  wind  speed. 
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Particle  diameter  (microns) 

Hgure  4-1.  Surface  deposition  and  gravitational  settling  velocity  versus  particle 
diameter  for  two  surface  roughness  lengths.  Surface  Motion  velocity 
is  assumed  to  be  0.4m/s,  and  the  deposition  is  calculated  for  a  height 
of  IQm. 


4JL3  Initial  Conditions. 

The  numerical  integration  is  initialized  with  the  dust  distribution  from  a  DICE 
calculation  at  f  =  5  min.  for  a  IMT  burst  at  a  scale  height  of  200sfL  The  DICE 
computation  (Hassig  et  al.,  1992,  p25)  was  made  with  a  very  fine  grid  near  the  surface,  in 
order  to  resolve  the  sweep-up  of  the  pedestal  dust  and  follow  its  evolution  for  the  first  S 
minutes.  The  DICE  computation  was  also  concerned  with  the  lar^r  size  ejecta,  so  that 
only  one  particle  size  group  was  devoted  to  the  sub-lOOpm  sizes,  and  at  5  minutes  this  is 
the  only  group  with  significant  mass  remaining  in  the  air.  The  spatial  distribution  of  the 
mass  in  this  group  is  given  in  Figure  4-2,  which  shows  the  cloud  extmiding  up  to  about 
15km,  with  a  shallow  pedestal  out  to  about  Skm  horizontally.  The  primary  focus  of  the 
DICE  calculation  was  to  examine  the  near-surface  dust  layer,  and  the  increased  diameter 
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Height 

(km) 


Figure  4-3.  Vertical  profiles  of  horizontally-integrated  dust  concratration  from 
the  DICE  run  at  tsSmins.  (a)  ^ows  full  profile  (note  logarithmic 
concentration  scale),  (b)  shows  detail  in  ^  lowest  levels. 
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of  the  cloud  stem  above  the  50Qm  height  is  due  to  increased  diffusion  -  both  numerical 
and  parameterized  turbulence.  Above  about  200m,  both  vertical  cell  size  and  turbulence 
mixing  length  increase  with  altitude.  The  horizontally-integrated  mass  distribution,  used 
to  initialize  the  boundary  layer  calculation,  is  givmi  in  Figure  4-3.  It  can  be  seen  that  the 
highest  densities  are  near  the  surface,  but  the  extent  of  the  high  density  is  only  about 
100m  in  the  vertical.  The  dust  distribution  does  extend  high  into  the  atmosphere, 
however,  and  over  the  long  period  of  the  integration  this  dust  can  influence 
concentrations  in  the  boundary  laym-  through  ^vitational  settling. 

We  avoid  specifying  detailed  size  distributions  at  this  stage  by  using  a  normalized 
dust  density  for  a  range  of  particle  size  bins.  In  fact,  we  initialize  7  bins  spanning  the 
range  from  litm  to  100pm  with  logarithmic  spacing.  The  particle  size  bins  evolve 
independently,  i.e.  the  concentration  in  one  bin  is  not  affected  by  concentrations  in  any 
other  bin.  since  we  ignore  any  change  in  particle  size  due  to  agglomeration  or  fracture. 
We  can  therefore  renormalize  the  dust  concentrations  to  simulate  any  desired  initial  size 
distribution,  subject  only  to  the  constraint  that  the  total  initial  mass  matches  the  DICE 
mass  field. 

The  most  commonly  used  size  distribution  is  the  a"*  spectrum,  where  a  refers  to 
the  particle  radius.  This  power  law  applies  to  the  number  density,  so  that  the  actual  mass 
spectrum  follows  an  cr^  law,  smce  the  particle  volume  is  proportional  to  a  ^ .  Most  of  our 
results  will  be  based  on  this  distribution,  but  we  will  examine  the  importance  of  the 
particle  size,  providing  guidance  on  the  sensitivities. 

4,2.4  Meteorological  Parameters. 

Results  were  obtained  for  a  range  of  meteorological  conditions,  in  an  effort  to 
gain  some  information  on  the  possible  range  of  variations  in  the  dust  field.  Ambient 
conditions,  such  as  wind  speed,  surface  heating,  and  surface  roughness  are  expected  to 
play  a  role  in  the  dust  field  evolution,  as  is  the  time  of  release  in  relation  to  the  diurnal 
variation.  We  have  computed  a  matrix  of  cases;  two  geostrophic  winds  were  considered, 
lOms'^  and  2ms*^  and  two  surface  roughness  lengths  were  used,  a  Im  roughness 
representative  of  a  forest-type  surface,  and  a  1cm  roughness  typical  of  a  anooth  desert- 
type  surface.  Two  surface  heating  rates  were  also  considered,  a  high  value  of 
0.275®Kms“i  (approximately  SOOWm"^)  and  a  low  value  of  0.075°Kms''*;  these  ate 
maximum  daytime  heating  rates,  with  the  high  value  bmng  typical  of  a  clear  summer  day. 
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and  the  low  value  more  appropriate  for  a  cloud-covered  situation.  For  all  the  cases 
considered,  four  initializations  were  made  at  6  hour  intervals  from  midnight  to  determine 
the  effect  of  time  of  day  of  release  on  the  dust  evolution. 

The  surface  heating  is  specified  as  a  sinusoidal  variation  during  the  daylight 
hours,  Le.  from  06:00  to  18:00  local  time,  with  a  small  constant  negative  heat  flux  at 
night  We  have  not  considered  variations  in  the  length  of  the  day,  since  we  expect  the 
specified  parameter  variations  to  dominate  the  response.  The  periodic  surface  heating 
allows  the  calculation  to  be  run  for  several  days,  so  that  a  near-periodic  diurnal  boundary 
layer  can  be  developed.  The  boundary  layer  is  not  perfectly  periodic  because  the  net  heat 
input  during  a  24  hour  period  is  positive,  te  more  heat  is  input  during  the  day  than  is  lost 
at  night  This  imbalance  causes  the  mixed  layer  to  grow  a  little  deeper  each  day  as  it 
steadily  warms  and  erodes  the  overlying  stable  atmosphere.  The  relative  change  from 
one  day  to  the  next  is  small,  however,  and  does  not  significantly  affect  the  dust  transport 
The  actual  surface  energy  balance  at  night  is  a  complex  interaction  between  turbulent 
flux,  radiative  fluxes  and  heat  flux  from  the  soil  or  subsurface,  which  determines  the 
nocturnal  boundary  layer  structure.  We  have  simply  Reified  a  constant  cooling  at  night 
however,  but  the  variation  of  external  parameters  considered  covers  the  range  of 
nocturnal  boundary  layer  behavior. 

All  the  boundary  layer  integrations  were  initialtzed  with  an  idealized  wind  profile, 
but  were  run  for  3-5  days  to  allow  the  establishment  of  a  reasonably  periodic  diurnal 
profile.  The  dust  releases  were  made  into  the  established  boundary  layer  at  00:00, 06:00, 
12:00  and  18:00  local  time,  where  12:00  represents  maximum  surface  heat  input  and 
06:00  and  18:00  represent  sunrise  and  sunset  respectively.  A  typical  variation  of  the 
vertical  turbulence  energy  component  for  a  2  day  period  is  shown  in  Figure  4-4.  These 
results  were  obtained  following  a  5  day  initialization  period,  and  clearly  show  the 
suppression  of  turbulence  during  the  night  Shortly  after  sunrise,  when  the  surface  heat 
flux  becomes  positive,  a  relatively  slow  growth  of  the  turbulence,  both  in  depth  and 
intensity  is  visible,  as  the  heating  erodes  the  surface-based  inversion  formed  during  the 
night  At  about  10:30  local  time,  the  mixing  reaches  the  top  of  the  nocturnal  inversion 
and  grows  very  quickly  through  the  remnant  of  the  previous  day's  mixed  layer  to  reach  a 
depth  of  about  1200m.  The  boundary  layer  grows  very  slowly  for  the  rest  of  the  day, 
with  slowly  decaying  intensity  after  12:00  until  the  surface  heating  becomes  negative 
about  18KX)  and  the  turbulence  colliq)ses.  This  behavior  is  closely  repeated  on  the  second 
day  with  a  boundary  layer  roughly  100m  deeper. 
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Hgure  4-4.  Diurnal  variation  of  the  planetary  boundary  layer. 

(a)  r.m.s.  vertical  velocity  fluctuation,  contour  interval  O.lm/s 

(b)  surface  heat  flux  in  KnJs 


All  the  mixing  la]«rs  used  in  diis  study  were  between  1000m  and  1500m  deep,  which  is  a 
typical  range  for  convectively  mixed  layers.  It  is  certainly  possible  for  deeper  mixing  to 
occur,  but  2000m  is  a  reasonable  extreme  value  for  inversion-c^q)ed  mixing.  For 
shallower  layers,  it  is  possible  for  large-scale  subsidence  to  prevent  mixed  layer  growth. 
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but  layeis  below  300m  are  very  unusual  under  conditicms  of  positive  surface  heating. 
These  extreme  values  would  modify  die  dust  concentration  results  somewhat,  so  that  a 
well-mixed  value  would  be  proportionally  higher  in  a  shallow  layer. 


43  RESULTS. 

Hgure  4-5  shows  the  concmitration  profiles  at  4, 8,  and  12  hours  after  release  for 
the  four  initial  start  times  with  Uo  -  10ms~^,  zo  ==  0.01m,  and  =  0.07S**Kms~^ 
("•100Wm~2).  The  initial  particle  size  distribution  is  cr^ ,  and  the  multiple  lines  indicate 
cumulative  total  up  to  the  various  size  groups,  which  are  represented  by  nominal  particle 
diameters  of  1,  2,  6,  10,  20,  60,  100^.  These  particle  sizes  can  be  considered  as  the 
masR  mean  diameters  (MMD)  of  7  size  bins,  with  boundaries  at  1, 1.4,  3.5, 7.8,  14,  35, 
78,  and  144pm.  The  MMD  of  the  smallest  size  group  is  actually  1.2pm.  not  1pm,  but  the 
size  distribution  is  truncated  at  1pm  and  the  discrepancy  is  not  significant  for  the 
deposition  processes  considered  in  this  sectioiL  Thus  the  leftmost  line  in  Hgure  4-5  gives 
concentration  in  the  1pm  group  (between  1  and  1.4pm),  the  next  is  concentration  in  the 
1pm  and  2pm  groups  (between  1  and  3.5pffl),  while  the  rightmost  line  gives  total 
concentration  for  all  size  groups.  The  boundary  layer  mixing  during  the  daytime  hours  is 
very  evident  in  this  figure,  which  shows  a  well  mixed  layer  reaching  depths  of  about 
1.5km  in  the  afternoon.  At  4  hours  after  release,  three  of  the  four  cases  show  roughly  the 
same  maximum  level,  namely  80-100  Mg/m,  but  it  is  clear  that  the  noon  release  has 
produced  a  deep  mixed  profile  and  reduced  the  mean  concentration  by  about  a  factor  of  2. 
The  other  release  times  are  confined  closer  to  the  surface  and  the  initial  concentration  is 
only  modified  by  the  gravitational  settling  of  the  larger  particles. 

At  later  times  the  boundary  layer  mixing  has  made  its  appearance  in  more  of  the 
cases,  so  that  at  12  hours  after  release  there  is  rfrong  similarity  betwemi  the  00:00, 06:00 
and  12:00  release  cases.  The  18:00  release  evolves  through  the  nocturnal  hours  during 
this  period,  and  is  still  shows  a  shallow  maximum,  although  the  turbulence  induced  by 
the  lOms*^  wind  does  mix  the  material  up  to  about  300m  depth.  The  other  three  cases 
have  all  undergone  at  least  one  half  of  the  daily  mixing  cycle,  and  show  very  similar 
profiles  with  an  almost  constant  total  concentration  of  about  35Mg/m.  The  12:00  release 
does  show  a  more  pronounced  surface  minimum  at  12  hours  since  the  deposition  removes 
the  lowest  layers  of  dost  and  there  is  no  replenishment  from  above  due  to  the  limited 
vertical  extmit  of  the  turbulence  in  the  nocturnal  hours.  The  preferential  loss  of  the  larger 
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(a)  Release  time  00:00 


(b)  Release  time  06:00 


(c)  Release  time  12:00 


(d)  Release  time  18:00 


Hgure  4-5.  Integrated  dust  profiles  at  4  8,  and  12  hours  after  release.  Release  is  at  (a)  00:00; 
(b)  ^:00;  (c)  12.-00;  (d)  18:00.  Wind  speed  lOm/s,  roughness  length  O.Olm 
(Continued). 
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pirticle  groups  is  also  evident  in  the  results  at  12  hours,  with  the  firaction  of  the  dust  in 
the  rightmost  bands  being  reduced.  It  is  clear,  however,  that  most  of  the  reduction  occurs 
in  the  first  4  hours,  as  the  larger  particles  fall  out  relatively  quickly.  Almost  all  of  the 
mass  is  contained  in  the  lowest  5  size  groups,  i.e.,  panicle  diameters  less  than  3Spm,  after 
4  hours  of  evolution. 

Increasing  the  surface  roughness  to  Im.  but  holding  all  other  factors  tlv;  same, 
produces  the  dust  profiles  in  Figure  4-6.  There  are  two  main  differences  induced  by  this 
change;  fiirst  the  nocturnal  boundary  layer  is  much  deeper,  so  that  the  pedestal  is  mixed 
higher  at  all  times,  and  second  the  deposition  is  larger,  so  that  the  concentrations  in  the 
lowest  1km  at  12  hours  after  release  are  slightly  less  than  those  in  Figure  4-S.  The 
differences  are  not  large,  however,  since  the  bulk  of  the  remaining  material  is  in  the  small 
size  groups,  which  stiU  have  small  deposition  velocities.  The  surface  minimum  is  more 
pronounced  in  this  case,  indicating  the  stronger  removal  mechanism  at  the  rough  surface. 
The  nocturnal  turbulence  is  generated  by  the  relatively  high  wind  speeds  in  Figure  4-6. 
This  can  be  seen  in  Figure  4-7,  which  shows  the  same  high  roughness  length,  but  has  a 
wind  speed  of  only  2ms~^  In  this  low  wind  case,  the  turbulence  is  quenched  almost 
down  to  the  surface  at  night,  so  that  the  initial  profile  evolves  over  12  hours  due  only  to 
gravitational  settling  in  the  18:00  release  case.  For  other  release  times,  however,  the 
daytime  mixing  produces  profiles  at  12  hours  which  are  similar  to  the  high  wind  speed 
case,  but  with  slightly  higher  values  due  to  the  reduced  turbulent  deposition  rate. 

Some  indication  of  the  range  of  profiles  is  given  in  Figure  4-8,  which  shows 
profiles  at  12  hours  after  an  06:00  release  under  various  meteorological  conditions.  All 
the  cases  show  a  mixed  layer  development,  with  low  level  concentrations  ranging  from 
roughly  20  Mg/m  to  35  Mg/m.  This  narrow  range  is  representative  of  the  complete 
extent  of  the  results  from  our  matrix  of  release  times,  wind  speeds,  etc.  The  high  wind 
speed,  high  roughness  length  case  shows  visible  evidence  of  surface  deposition,  with  a 
detectable  low  level  minimum.  The  relatively  small  variation  with  meteorological 
parameters  is  due  to  the  similarity  of  the  mixing  history  over  a  12  hour  period  between 
the  various  cases,  together  widi  the  small  deposition  rates  for  tte  small  particles.  The 
general  conservation  of  the  small  particle  mass  ensures  that  the  mixed  layer  density  is 
always  about  the  same,  and  since  all  our  mixing  layos  are  roughly  the  same  depth  we 
obtain  the  same  final  concentration. 


41 


(a)  Release  time  00:00 


(b)  Release  time  06:00 
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Hgure  4-6.  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  (a)  00:0 
0>)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  lOm/s,  rougtoess  length  1.0m. 
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(c)  Release  time  12:00 


(d)  Release  time  18:00 
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Figure  4-6.  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  (a)  00:00; 
(b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  lOm/s,  roughness  length  l.Cto 
(Continu^). 
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Figure  4-7.  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  (a)  00:00; 
(b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  2in/s,  roughness  length  1.0m. 
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Figure  4-7  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  U)  00:00; 
(b)  06:00;  (c)  12:00;  (d)  18:00.  Wind  speed  2in/s,  roughness  length  1.0m 
(Continu^). 


The  contribution  from  the  pedestal  in  the  late  time  dust  profiles  is  indicated  in 
Figure  4-9,  where  the  initial  dust  field  was  truncated  at  z^sOOm  but  all  other  parameters 
wme  the  same  as  those  in  Figure  4-5.  The  limited  vertical  extent  is  evidmit  in  Figure  4-9, 
with  no  dust  present  above  die  mixed  layer,.  The  levels  are  also  reduced  at  12  hours, 
with  mixed  layer  levels  of  about  20  Mg/m  in  the  mixed  cases  and  about  60  Mg/m  in  the 
IS.KX)  release  case,  Le.  roughly  2/3  of  the  full  vertical  distribution  results.  We  also  note 
that  the  largest  size  groups,  with  MMD  60pm  and  100pm,  are  completely  removed  from 
die  pedestal  by  12  hours,  and  that  the  100pm  group  has  settled  out  before  4  hours.  The 
role  of  the  tuibulence  in  lofting  the  larger  particles  can  be  seen  in  the  presence  of  the 
larger  particles  in  the  upper  region  of  the  mixed  layer  4  hours  after  the  12:00  release, 
which  immediately  experiences  turbulent  mixing,  but  not  for  any  of  the  other  release 
times. 


The  particle  size  is  very  important  in  determining  the  fate  of  the  dust  over  periods 
of  12  hours  or  more.  Figure  4-10  shows  the  same  cases  as  Hgure  4-8,  but  under  the 
assumption  that  all  the  particles  are  10pm  size,  instead  of  the  spectral  distribution.  For 
high  wind  speed  cases,  the  concentrations  at  12  hours  are  lower,  indicating  increased 
deposition  and  gravitational  settling  compared  with  the  earlier  figure.  The  high 
roughness  and  high  wind  speed  results  are  interesting,  since  they  demonstrate  that 
turbulent  dqx>sition  at  the  surface  is  the  important  effect  in  reducing  the  low  level  dust 
concentrations.  Lx)w  wind  speed  cases  show  slightly  higher  concentrations  than  Figure 
4-8.  This  could  be  anticipated  from  Figure  4-1,  which  suggests  that  the  deposition  rate 
for  these  particles  is  important  in  the  high  wind  speed,  high  roughness  cases..  The 
deposition  velocity  is  generally  small  for  1pm  particles,  but  10pm  particles  can  be 
absorbed  quite  rapidly  under  the  rq)proptiate  conditions  of  high  turbulence  and  roughness. 
The  reduction  in  concmitration  near  the  surface  is  also  evident  in  the  high  wind  speed 
cases,  indicating  the  importance  of  the  surface  deposition.  In  general,  these  results  show 
much  mote  sensitiviQr  to  the  meteorology  than  tiiose  in  Figure  4-8. 

A  mote  realistic  particle  size  distribution,  the  proposed  incipient  distribution  for 
tiie  sweep-up  material  (Hassig,  personal  communication),  emphasizes  the  mass  in  the 
larger  particles.  Results  for  various  meteorological  conditions  are  given  in  Figure  4-11. 
There  is  clearly  more  material  in  the  20pm  and  60pm  size  groups,  and  this  material 
di^lays  the  expected  sensitiviQr  to  meteorology  indicated  by  Hgure  4-10.  In  general,  the 
increased  proportion  of  mass  in  the  larger  size  groups  leads  to  lower  concentrations  at 
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Figure  4-9.  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  (a)  00:00; 
(b)  06:00;  (c)  12KX);  (d)  18:00.  Wind  speed  Khn/s,  roughness  lengdi  0.01m, 
pedestal  dust  only. 
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Figure  4-9.  Integrated  dust  profiles  at  4, 8,  and  12  hours  after  release.  Release  is  at  (a)  00:00; 
(b)  06KX);  (c)  12:00;  (d)  18:00.  Wind  q;)eed  Khn/s,  roughness  length  0.01m, 
pedestal  dust  only  (Continued). 
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Figuie  4-10.  Integrated  dust  profiles  at  12  hours  after  06:00  release  assuming  all  dust  is 
10  microns.  Surface  heat  flux  maximum  is  0.07SKm/s  unless  noted. 
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Hgure  4-11.  Integrated  dust  profiles  at  12  hours  after  06:00  release  assuming  proposed  incipient 
particle  size  distribution.  Surface  heat  flux  maximum  is  0.075KWs  unless  noted. 
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later  tunes,  since  this  mass  is  dq>osited  more  effectively  than  the  smallest  size  groups. 
The  mixed  layer  total  concentrations  are  roughly  half  die  values  of  Hgure  4-8. 

4.4  SUMMARY  OF  BOUNDARY  LAYER  RESULTS. 

Numerical  solutions  of  the  turbulent  transport  and  deposition  model  equatitHis 
have  indicated  that  a  significant  fraction  of  the  dust  pedestal  from  a  large  nuclear  burst 
can  remain  lofted  in  the  surface  mixing  layer  for  periods  of  up  to  24  hours.  This  behavior 
is  dependent  on  the  particle  size  distribution  in  the  pedestal,  and  can  be  summarized  as 
follows. 

For  particles  below  Spm,  the  turbulent  deposition  and  also  gravitational  settling 
are  msgligible  over  the  time  period  of  one  day,  so  that  the  initial  material  is  simply 
redistributed  vertically  by  the  turbulence  induce  by  daytime  heating.  Provided  a  mixed 
layer  of  depth  around  1km  is  formed  during  the  day,  the  vertical  distribution  of  material 
is  insensitive  to  meteorological  parameters,  although  the  lateral  distribution  will  depend 
on  the  wind.  This  result  will  apply  if  the  bulk  of  the  pedestal  dust  is  in  this  size  range,  as 
might  be  predicted  by  a  typical  size  spectrum  assumption,  if  a  lower  cutoff  of  1pm  is 
assumed. 

For  particles  above  30pm,  the  gravitational  settling  is  sufficiently  rapid  that  the 
pedestal  is  removed  after  a  few  hours,  even  in  the  preseix:e  of  relatively  strong  vertical 
nuxing.  This  result  is  also  insensitive  to  meteorological  variations. 

The  size  range  between  5pm  and  30pm  di^lays  more  complex  behavior,  since 
the  turbulent  deposition  velocities  in  this  range  can  vary  from  negligible  to  significant 
over  the  range  of  meteorological  parameters,  while  the  gravitational  settling  is  generally 
smalL  Thus,  under  light  wind,  small  roughness  conditions  most  of  this  material  will 
remain  lofted  through  the  mixed  h  for  24  hours,  while  the  surface  deposition  under 
high  wind  speed,  high  roughness  col  litions  will  deplete  the  material  by  up  to  90%  over 
such  a  period. 
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SECTION  5 


STATISTICAL  FALL-OUT  ANALYSIS 

The  dust  lofted  by  a  nuclear  burst  persists  in  the  atmosphere  for  many  hours. 
Different  size  particles  fall  at  different  rates,  so  the  distribution  changes  with  time  as  it  is 
moved  by  the  local  winds.  Dust  can  cause  damage  to  aircraft  systems,  so  [wedictions  of 
the  late-time  distribution  are  important  for  planning  purposes.  However,  while 
predictions  ate  certainly  feasible  for  a  given  wind  field  using  current  late-time  dispa:sion 
models,  the  wind  field  is  genoally  not  known  for  planning  scenarios.  Statistical  studies 
of  the  climatological  variability  in  the  wind  field  at  specific  sites  shows  the  r.m.s. 
fluctuation  in  wind  component  speed  to  be  larger  than  the  mean  value,  indicating  a  large 
uncertainty,  and  a  statistical  description  of  the  effects  of  this  variability  on  aircraft 
vulnerability  is  desirable. 

R(^  (1991)  has  initiated  an  effort  to  extend  the  SAFER  methodology  (Ross  and 
Mazzola,  1993)  to  the  aircraft  problem.  This  approach,  successfully  applied  to  the  RV- 
fratiicide  problem,  gives  a  very  efficient  analytical  description  of  the  probabili^  of 
encountering  critical  loading,  and  is  simple  enough  to  iqiply  at  a  systems  planning  level 
The  difficulties  involved  in  the  late-time  aircraft  problem,  however,  prevented  the 
calculation  of  more  than  the  ensemble  average  dost  mass  loading  along  a  fligh^)ath.  The 
late-dme  problem  is  more  complicated  than  the  RV-fratridde  problem,  since  it  critically 
involves  dust  falling  through  different  levels  of  the  atmosphere,  and  the  winds  vary  vdth 
altitude  in  general  The  RV-fcattidde  problem  was  restricted  to  relativdy  short  times, 
where  only  horizontal  transport  need  be  considered.  Unfortunately,  the  ensemble  mean 
mass  loading  is  not  the  most  useful  quantity  ftn*  a  planner,  who  usually  wants  to  know  the 
probability  of  failure,  and  the  objective  of  the  present  study  is  the  extension  of  Ross' 
tqpproach  to  provide  this  probabilistic  information. 

The  dust  density  distribution  is  specified  as  a  function  of  space  (2).  time  (r),  and 
particle  size  (a),  le.,;ir(x,r,a).  The  aircraft  mass  loading  involves  a  line  integral  of  the 
dust  along  a  flightpath  at  a  certain  altitude.  This  could  be  in  any  direction  at  any  time 
after  the  burst,  but  we  specify  tlmse  parameters  as  fixed  for  definiteness  in  the  analysis. 
We  considCT  the  mass  loading  along  a  N-S  flightpath  at  z  =  100m  and  t  =8  hours  after 
burst. 
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(5.1) 


m(x)  5s  Jrfy  “Yda  ;f(ac,y,100in,  8hrs,  a) 
=  ‘J^l{x,a)da 


where  /u(x,a)  is  the  mass  loading  density  from  particles  size  a.  We  assume  that  nixxi)  is 
obtained  from  the  initial  loading  of  particle  size  a.  at  an  initial  altitude  of  zo(a)  such  that 
tlm  particles  fall  to  z  s  100m  after  8  hours.  This  is  the  methodology  adopted  by  Ross 
(1991)  and  provides  a  simple  conceptual  frameworit  for  analyzing  the  effects  of  wind 


variations. 


We  further  assume  that  /i(x^)  has  a  Gaussian  spatial  distribution  with  centroid 
and  spread  (3(a),  i.e., 


,  .  M^(a)  r  (x-X(a)f‘ 


where 


M>(«)  =  llz(x*y,Zo(a),0,a)  dxtfy 


The  centroid  location,  X,  is  determined  by  the  wind  profile,  u{z),  over  the  levels 
travmsed  by  the  particles  as  they  fall  to  z  »  100m.  Assuming  a  constant  fall  velocity, 

X  =  (/(zo)r  (5.4) 


where 


(z~100m)„i„ 


(5.5a) 


l.e.,  the  centroid  displacement  is  the  layer-averaged  velocity  multiplied  by  the  translation 
time.  In  reality,  the  particle  fall  speed  is  height-dependent  due  to  the  variation  in  air 
density,  so  the  layer-average  velocity  should  strictly  include  a  weighting  factor 
proportional  to  the  inverse  fall  speed,  as  described  by  Ross  (1991).  Thus 


l/(z)  =  — 

8hrsjw^(z0 


(5.5b) 


where  is  the  local  fall  ^>eed.  Definition  (5.5b)  is  size-dependent,  in  general,  since  fall 
speed  depends  on  particle  size,  so  the  wind  statistics  will  be  defined  using  Equation 
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(5  Ja)  in  the  subsequent  analysis  for  simplicity.  We  note  that  the  proper  definition  (S.5b) 
could  be  used  to  define  the  appropriate  l/(z),  and  this  was  actually  used  in  the  Monte- 
Carlo  realizations  computed  below. 

The  spread,  <7,  is  determined  by  atmospheric  diffusion.  We  assume  a  simple 
model  for  diffusion,  although  a  more  sophisticated  formulation  could  easily  be  employed 
since  the  subsequent  analysis  only  uses  the  value  of  a  and  is  independent  of  the  precise 
diffusion  mechanism.  For  definiteness,  we  use  the  Heffter  parameterization  for  the 
diffusive  spread,  which  is  a  constant  rate  growth  of  the  horizontal  spread.  i.e., 

<s  =  0.5r  ((T  in  meters,  t  in  seconds)  (5.6) 

Equations  (S.l)-(5.6)  permit  the  calculation  of  m(r)  if  we  are  given: 

(i)  a  wind  profile,  and 

(ii)  an  initial  dust  (tensity  distribution. 

The  problem  we  seek  to  solve  is  to  determiite  the  statistics  of  m(x)  when  the  winds  are 
uncertain,  i.e.,  we  only  have  statistical  information  on  the  wind  field.  However, 
Equations  (5.1)-(S.6)  can  be  used  directly  to  perform  Monte-Carlo  simulations  of  the 
loading  problem  using  a  large  collection  of  observed  wind  profiles.  This  technique  is  not 
computationally  efficient,  but  does  generate  a  'true'  distribution  with  which  to  compare 
simplified  models. 

The  main  component  of  uncertainty  in  m  is  due  to  uncertainty  in  X,  i.e.,  where 
does  the  dust  go?  We  therefore  assume  tris  fixed  in  the  following  analysis,  and  only  the 
centroid  is  random.  The  conceptual  model  is  thus  a  horizontal  Gaussian  spatial 
distribution  with  fixed  spread  for  each  particle  size  located  at  a  random  centroid  positioiL 

The  mass  centroid,  X,  is  a  random  variable,  related  to  the  random  wind  field 
through  Equation  (5.4).  The  mean  and  variance  of  X  are  easily  seen  to  be 

X  =  (7(zo)f  (5.7) 

r=jF  =  (X-X)' =i?^(z,)f^  (5.8) 

We  represent  the  probability  distribution  function  (pdf)  of  X  by  a  Gaussian,  so  that 


Then,  using  Equation  (5.2),  we  calculate  the  mean  size-dependent  loading. 
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]i{x,a)^]tiix,a;X)piX)dX 

(5.10) 

where  =  +  Z’ 

and  hence,  from  Equation  (5.1),  the  mean  mass  loading  is 

(5.11) 

mix)  =  *Jjiix,a)da 

(5.12) 

•l 


The  mean  wind  statistics  are  presented  by  Ross  (1991)  and  are  repeated  here  in 
Figure  5-1,  which  shows  the  mean  and  standard  deviation  as  a  function  of  altitude  from 
3652  Moscow  wind  profiles.  These  profiles  were  used  to  compute  3652  mass  loading 
profiles  from  an  initial  TASS  LM02  cloud  using  Equations  (5.1)-(5.6).  The  ensemble 
mean  loading  from  these  Monte-C^o  results  is  shown  in  Figure  5-2,  along  with  the 
ensemble  prediction  from  equations  (10)  and  (12).  The  mean  loading  compares  well  with 
the  Monte-Carlo  calculation,  indicating  the  accuracy  of  the  Gaussian  approximation  for 
the  pdf  of  the  wind  fluctuations. 

The  above  analysis  gives  qualitatively  amilar  results  to  Ross's  analysis,  in  spite  of 
his  assumption  of  disks  rather  than  Gaussian  distributions.  The  main  difference  is  due  to 
our  assumption  that  the  wind  is  constant  below  the  lowest  observation  level,  whereas 
Ross  tq)pears  to  assume  that  the  wind  goes  linearly  to  zero  at  the  ground.  We  believe  that 
the  mixing  effect  of  the  planetary  boundary  layer  will  produce  a  more  uniform  vertical 
profile  over  this  region.  We  also  include  the  late-time  diffusion  using  the  Heffter  model, 
so  our  intercepted  dust  loadings  are  lower  than  Ross'.  The  mean  mass  loading  is  not 
sensitive  to  the  detailed  shape  assumption,  since  a«Z  ,  and  therefore  the  wind 
fluctuations  dominate  the  result 
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Figure  5-1.  Wind  statistics  for  the  East-West  component  firom  3652  Moscow  wind  profiles, 
(a)  local  wind  at  altitude;  (b)  average  wind  between  altitude  and  ground. 
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Hgure  5-2.  Mean  mass  loading  after  8  hours  along  a  North-South  lligh^atfa  at  100m 
altitude  for  an  ensemble  of  3652  Moscow  winds.  The  explicit  Monte- 
Carlo  result  (solid)  and  the  ensemble  statistical  calculation  (dashed)  are 
shown  as  a  fimction  of  East-West  location  relative  to  the  burst 
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So  far,  we  have  only  produced  a  slighdy  different  version  of  the  mean  mass 
loading.  CXir  objective,  howev^,  is  to  obtain  information  about  the  pdf  of  mCx),  not  just 
the  average  value,  so  we  need  to  extend  die  analysis.  Calculation  of  the  detailed  pdf  of  m 
is  extremely  complicated,  so  a  simpler  i^proach  is  needed.  The  mean  value  is  the  first 
moment  of  the  pdf,  and  gives  no  information  about  the  range  of  possible  values  in  the 
distribution.  By  calculating  higher  moments,  more  detail  about  the  pdf  is  revealed  and 
quantitative  estimates  can  be  derived.  We  seek  to  compute  the  second  moment,  m'^ ,  and 
use  an  analytic  shape  assumption  about  the  pdf  to  get  estimates  of  likelihood  of 
exceeding  a  specified  threshold.  The  anticipated  sh^  is  the  clipped-normal  function, 
which  provides  an  intermittent  pdf  and  is  completely  determined  by  the  mean  and 
variance  of  the  distribution. 


Some  support  for  the  use  of  a  clipped  normal  distribution  is  provided  by  the 
Monte-Carlo  results  themselves.  The  explicit  calculation  of  3652  realizations  allows  us 
to  construct  the  pdf  at  any  location,  and  tte  cumulative  distribution  function  (cdf)  of  the 
mass  loading  at  two  x-locations  is  shown  in  Hgute  S-3.  The  cdf  represents  the  fraction  of 
the  loadings  below  the  level  plotted  on  the  horizontal  axis.  The  distortion  of  the  vertical 
scale  is  such  that  a  clipped  normal  distribution  is  a  straight  line,  and  the  intercept  on  the 
axis  gives  the  probability  of  zero  loading,  Le.  the  intermittency.  The  Monte-Cario  results 
show  a  reasonably  straight-line  section,  apart  from  a  tail-off  at  small  loadings,  i.e.,  there 
is  a  sigiuficant  probability  of  a  small  but  non-zero  loading.  The  pdf  at  higher  loadings 
can  be  well  represented  by  the  clipped  normal,  but  the  details  at  low  mass  loading  would 
require  a  different  shq)e  assumption  for  the  pdf.  For  the  present  purposes,  we  assume 
that  the  larger  values  are  of  interest,  and  use  the  clipped  normal  distribution. 


The  statistics  of  nix.a)  are  identical  to  those  of  Gifford's  (1959)  meandering 
plume  model  Thus,  following  Gifford’s  analysis,  the  variance  is  given  by 

^  =  j(M-JD^p(X)dX 


Ml 


Inaifp- 


exp 


ix-xf 

(t"+2I\ 


(5.13) 


However,  the  mass  loading  involves  an  integral  over  particle  size,  so 
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Figuie  5-3.  Cumulative  distribution  of  tl»  mass  loading  at  8  hours  from  the  Monte- 
Cario  results  with  3652  winds  (Le.  probabiUty  that  loading  is  less  than 
given  value)  at  two  locations. 
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A  similar  problem  arises  in  the  time-averaging  of  Gifford's  meandering  plume, 
where  temporal  correlations  ate  required,  and  a  solution  based  on  the  wind  fluctuation 
correlation  was  given  by  Sykes  (1984).  We  can  use  the  same  ^proach  for  our  particle 
size  integration,  but  we  need  to  extend  the  analysis  to  account  for  the  inhomogeneous 
statistics  of  the  wind  field.  As  we  shall  see  below,  the  wind  field  is  very  strongly 
correlated  at  different  heights,  but  the  increasing  mean  and  variance  with  altitude  has  a 
large  effect  on  the  loading  correlations.  The  two-point  correlation  is  defined  as 

/i(a,)/i(aa) «  (5.15) 

and  this  can  be  evaluated  if  the  joint  probability  of  the  centroid  locations,  for 

the  two  particle  sizes  can  be  ^lecified. 


Since  the  correlation  r  is  related  to  the  wind  correlation.  We  model 

r(fl,,flj)  =  exp 

wtere 


(5.17) 

(5.18) 
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and  u  determined  Crmn  dm  coneladon  of  the  winds  at  different  altitudes.  The 
coneladon  functicMu 


^  (5.19) 

is  shown  in  Hgure  S-4  for  ztsl450m.  The  very  strong  correlation  over  a  depth  of  15km 
is  obvious.  The  actual  velocity  correlation  is  also  given,  illustrating  the  effea  of  the 
vertical  averaging  in  (5-Sa).  The  fact  that  a  falling  particle  experiences  an  average 
velocity  over  the  depth  of  its  fall  gready  increases  the  correlation  between  the  winds 
e3q)erienoed  by  particles  fiom  different  initial  heights.  Htting  an  exptmential  fumction  to 
the  correlation  curve  in  Figure  S<4  gives  a  vertical  length  scale  of  about  Sfflon  for  the 
velocity  field.  Thus  the  transport  at  altitude  is  extremely  well  correlated  with  the  near- 
surface  transport,  but  it  is  not  identical  since  the  winds  are  higher  aloft  This  difference 
causes  different  particle  sizes  to  fall  out  at  different  locations  and  destroys  the  particle 
correlations. 


Substituting  (5-16)  into  (5-15)  and  using  (5-2),  after  some  algebraic  manipulation, 
we  obtain 


/x(ai)/i(a2)  =  Qb  exp 


(520) 


where 


Qa  ” 


_ ^01^02 _ 

2(^  +  ^-2ra^<72) 

(5.21) 
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Figure  S-4.  Correlation  between  East-West  velocities  at  two  altitudes  finom  3652 
Moscow  wind  profiles.  Lower  level  is  fixed  at  14S0bi.  Correlations 
between  the  lo^  velocity  and  the  average  over  the  layer  between  the 
two  levels  are  shown. 
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(5.22) 


(5.23) 


Equation  (S.14)  involves  a  double  integral  over  the  size  spectrum,  but  much  of 
this  space  is  uncorrelated.  For  speed  and  simplicity  of  computation,  we  represent  the 
correlation  function  by  an  exponential  with  scale  X,  and  integrate  analytically  over  one  of 
the  dimensions,  to  get 

m'^ix)  =  *7  ^'^(a,x)g(a,x)da  (5.24) 


udiete 


g(fl.x)  =  A(a.x)(2  -  + «-<--->'*) 


(5.25) 


The  correlation  scale,  A(a,x),  is  numerically  determined  from  the  point,  a',  where 

'  11,  I  ^  falls  to  a  value  of  0.5.  This  involves  a  numerical  search  at  each  point,  but 

H'Ha) 

is  sdU  very  efficient 


The  profile  of  the  mass  loading  standard  deviation  is  shown  in  Figure  5-5,  in 
comparison  with  the  Monte-Carlo  result  The  analytic  model  shows  reasonably  good 
agreement  with  the  explicit  calculatioit  giving  a  standard  deviation  roughly  twice  the 
mean  value  near  the  peak  of  the  curve  at  xsO.  As  we  move  away  from  the  peak,  die  ratio 
gets  larger  indicating  more  intermittent  statistics  as  illustrated  in  Hgure  5-3. 


Figures  S-6  and  S-7  show  the  predicted  probability  of  exceeding  Img/cm^  and 
lOmg/cm^  as  a  function  of  horizontal  (E-W)  location,  using  the  predicted  variance  and 
the  clipped  normal  assumption  or  the  pdf.  Again,  the  agreement  with  the  Monte-Carlo 
calculation  is  very  encouraging,  demonstrating  that  the  analysis  can  be  used 
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to  give  useful  statistical  predictions  of  the  probability  of  intercepting  various  mass 
loadings  in  the  presence  of  large  wind  uncertainties. 

The  details  of  the  statistics  can  be  easily  modified  to  accommodate  the  height- 
dependent  fall  speed  or  alternative  cloud  diffusion  rates.  The  analysis  could  also  be 
performed  for  subsets  of  the  climatological  statistics,  e.g.  seasonal  ensembles.  The 
representation  of  multiple  bursts  and  multiple  aircraft  routes  would  require  further  worit. 
but  the  method  presented  here  is  a  first  step  toward  a  realistic  assessment  of  the 
uncertainties  induced  by  the  wind  variability. 
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SECTION  6 


SCIPUFF  MODEL  STUDIES 

6.1  RADIOACTIVE  FALLOUT  FROM  20KT  WEAPON. 

A  calculation  of  the  late-dme  radioactive  cloud  and  fallout  from  a  small  nuclear 
weapon  was  performed  using  initial  conditions  from  the  DICE/MAZ  code.  The  standard 
cloud  description  was  extended  to  track  a  level  of  radioactivity  for  each  puff,  specified  as 
an  active  mass  fraction,/jt.  This  fraction  is  conserved  for  an  individual  puff,  assuming  no 
radioactive  decay  on  the  time  scales  of  interest,  and  the  total  active  mass  is  conserved  in  a 
puff  merging  operation.  In  addition  to  the  masr.  activity  level,  a  fallout  computation  was 
included  in  SCIPUFF  to  track  the  deposition  of  material  and  the  radioactivity  on  a  surface 
grid.  The  total  mass  fallout  from  each  puff  was  already  routinely  computed  at  each 
timestep,  since  this  is  the  mass  rate-of-change  term,  but  the  deposition  was  not  stored. 
The  code  was  extended  to  partition  the  mass  deposition  from  each  puff  on  a  surface  grid 
of  points,  using  a  Oaussian  distribution  in  the  horizontal  plane  and  an  appropriate  mass- 
conserving  normalization.  Thus,  all  mass  lost  from  the  puffs  is  balanced  by  an  increase 
in  the  surface  deposition. 

Initial  conditions  for  the  late-time  SCIPUFF  calculation  were  taken  from  the 
MAZ  results  for  a  2(NcT  surface  burst  at  r  =  5  minutes.  The  MAZ  calculation  tracked  the 
evolution  of  the  radioactive  isotopes  and  their  attachment  to  dust  particles,  and  was 
similar  to  the  demonstration  run  reported  in  Hassig  et  al.  (1992).  The  total  amount  of 
radioactive  material  released  in  the  explosion  is  ^proximately  1kg.  Hgure  6-1  shows  the 
vertical  distribution  of  the  radioactive  mass  for  particles  smaller  than  1mm,  and  also  the 
distribution  by  particle  size.  The  activity  is  confmed  to  the  main  cloud,  centered  at  6km 
altitude,  and  is  mostly  carried  by  the  particles  smaller  than  30|im.  There  is  a  small 
fraction  of  the  radioactivity  carried  by  particles  larger  than  1mm,  but  these  particles  ate 
deposited  very  close  to  the  burst  location.  The  subsequent  transport  and  diffusion  of  the 
smaller  particles  in  the  atmosphere  was  treated  as  a  passive  process,  with  the  exception  of 
gravitational  settling.  The  high  resolution  wind  fields  provided  by  the  ANATEX  data 
(see  Section  3)  were  used  to  transport  the  cloud  from  a  release  location  of  42*N  88*W. 
The  wind  vectors  at  four  altitudes  over  the  region  of  interest  are  displayed  in  Figure  6-2. 
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Figure  6-2.  Initial  wind  field  for  the  SCIPUFF  simulation 


The  winds  aloft  show  a  strong  Noitho'ly  flow  while  the  low  level  winds  near  the  release 
point  are  light  and  variable. 

The  mean  cloud  at  r  =  8  hours  and  two  altitudes  is  shown  in  Figure  6-3.  The  cross 
section  at  z  =  lOQm  shows  a  widely  spread  cloud  while  the  section  at  5km  shows  a  more 
compact  cloud  transported  further  South  over  Tennessee.  The  elongation  of  the  low 
altitude  cloud  is  a  result  of  the  strong  wind  shear  between  the  surface  and  Skm.  The 
standard  deviation  of  the  expected  concentration  value  is  also  indicated  in  Figure  6-3,  and 
is  comparable  with  the  mean  value.  The  accumulated  surface  dose  after  8  hours  in  p.gm~^ 
of  radioactive  material  is  shown  in  Hgure  6-4.  and  shows  the  same  elongation  toward  the 
South  as  the  mean  concentration  at  z=100m.  There  is  some  Northward  transport  in  the 
light  winds  around  the  source,  but  the  higher  altitude  winds  transport  the  larger  particles 
Southward  before  they  fall  to  the  ground. 

The  mean  radiation  dose  from  the  deposited  material  close  to  the  source  after  8 
hours  is  shown  in  Hgure  6-5.  We  have  converted  the  mass  deposition  to  equivalent 
radiation  at  Ihr  using  a  conversion  factor  of  10  rad/hr  per  l)igm~2  (R.  Christian,  personal 
communication).  The  high  doses  are  concentrated  in  a  plume  extending  to  the  South-East 
and  are  caused  by  the  fallout  of  the  larger  particles.  The  very  low  level  contamination  to 
the  North-West  in  Figure  6-3  results  from  small  particles  near  the  surface  with  little 
radioactive  materiaL  The  statistical  prediction  can  be  used  to  assess  the  probability  of 
exceeding  a  particular  dose  level.  This  uncertainty  is  a  result  of  the  incomplete 
knowledge  of  the  wind  field  as  represented  in  Figure  6-2  and  is  characterized  by  the 
concentration  fluctuation  variance.  Figure  6-6  shows  the  probability  of  exceeding  a  dose 
of  0.1  rad/hr  and  lOrad/hr  at  Ihr. 
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Hguie  6-3.  SOPUFF  prediction  of  the  total  dust  concentration  at  z=10Qni  and  z=Skm  at 
8  hours  after  initiation  of  a  2(8tT  bwst  centered  at  42^,  88°W. 

Contour  levels  of  10-3  10-2^  lO-i,  l  and  10  pg/nP. 
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Hgme  6-4.  SQPUFF  prediction  of  the  ensemble  mean  surface  dose  8  hours  afto- 
initiation  of  a  2(HcT  burst  centered  at  42^,88°W.  Contour  levels  of  lO'^, 
10^,  10■^  1(H,  lO-^and  IQ-^pg/n^, 


Figure  6-5.  SOPUFF  prediction  of  the  ensemble  mean  surface  radiation  dose  8  hours 

after  initiation  of  a  2(MtT  burst  Contour  levels  of  .01,  .1, 1, 10  and  100  rad/hr. 
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Figure  6-6.  SCIPUFF  prediction  of  the  probability  of  exceeding  a  surface  radiation  dose 
of  (a)  0.1  rad/hr  and  (b)  10  rad/hr  8  hours  aftor  initiation  of  a  2(&T  burst 
Contour  levels  of  1%,  5%,  20%,  50%.  80%  100%. 
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NUCLEAR  POWER  PLANT  RELEASE  SENSITIVITIES. 


In  the  event  of  a  nuclear  power  plant  explosion,  either  accidental  or  as  a  result  of 
intentional  attack,  the  release  of  radioactive  material  into  the  atmosphere  causes 
hazardous  conditions  for  human  activity.  The  hazards  range  from  intense  radiation 
exposure  close  to  the  source  all  the  way  to  large-area  and  long-term  health  effects  due  to 
radioactive  fallout  The  capability  to  predict  the  environmental  effects  of  a  nuclear  plant 
release  following  an  explosive  event  is  an  important  planning  tool,  involving  a  source 
estimate,  an  atmospheric  di^tersion  prediction,  and  a  health  effects  model  In  order  to 
assist  in  the  requirement  definition  for  the  source  calculation,  it  is  instructive  to  examine 
the  subsequent  atmospheric  processes  and  their  effects  on  the  radioactive  cloud.  The 
question  addressed  here  is  the  following;  how  detailed  a  description  of  the  initial,  local 
environment  do  we  need  to  be  able  to  reasonably  estimate  the  long  range  (lO-lOOOkm), 
long  term  (hours-days)  consequences?  SCIPUFF  has  been  run  for  a  variety  of  scenarios 
in  an  effort  to  assess  the  sensitivity  to  source  details.  The  matrix  of  SCIPUFF  runs 
included  variations  in  the  meteorology,  release  point,  and  source  characteristics.  The 
meteorology  variations  included  different  wind  fields,  release  points,  and  release  times. 

The  wind  fields  were  taken  from  tl^  National  Meteorological  Center’s  Nested 
Grid  Model  (NGM)  output  The  data  is  available  every  2  hours  for  the  first  quarter  of 
1987,  being  archived  as  part  of  ANATEX  (Across  North  America  Tracer  Experiment 
see  Section  3).  A  standard  boundary  layer  was  uniformly  imposed  over  the  entire  domain 
and  exhibited  a  diurnal  variation  with  a  depth  of  50m  during  the  night  and  a  maximum 
depth  of  1000m  during  the  day  in  January  and  February,  and  1500m  during  March. 
These  are  typical  values  for  the  planetary  boundary  layer  and  were  used  in  the  SCIPUFF 
validation  study  performed  with  the  ANATEX  data. 

Three  wind  field  data  sets,  each  covering  60  hours,  were  used: 

01/07/87:00:00  -  01/09/87:12:00 
02/04/87:00.00  -  02/06/87:12:00 
03/04/87:00:00  -  03/06/87:12:00. 

All  times  referred  to  in  this  study  are  C!entral  Standard  lime. 

Releases  were  made  from  two  locations,  central  Illinois  (90°W,  40°N)  and  the 
Four  Comers  area  of  Colorado  (109®W,  37®N)  for  both  a  midnight  (00:00CST)  and  noon 
(12:00CST)  explosion.  The  locations  were  chosen  to  give  a  broader  variation  of  wind 


77 


fields  and  the  release  times  give  a  smisidvity  to  diurnal  variations.  The  possible  range  of 
meteorological  conditions  is  obviously  very  wide  indeed,  but  this  limited  sample  should 
provide  some  insight  into  the  variability. 

The  puff  source  is  defined  by  five  parameters;  initial  cloud  size  and  height, 
duration  and  rate  of  rdease  and  release  material,  llie  size,  rate  and  duration  were  kept 
fixed  for  all  runs.  All  releases  were  for  12  hours  with  a  release  rate  of  1.0  kg/s, 
simulating  an  extended  release  at  constant  rate.  The  Gaussian  puffs  were  initialized  with 
a  spread  (o)  of  50m.  Four  release  heights  were  used:  200m,  SOOm,  100(hn,  and  2000ro. 
Four  different  source  materials  were  also  considmed  with  different  deposition  and  settling 
rates.  Settling  refers  to  the  equilibrium  fall  speed  under  the  influence  of  gravity,  while 
deposition  rate  is  the  speed  at  which  surface  absorption  occurs.  In  this  study,  we  have 
ignored  turbulent  deposition  of  particulates  but  used  a  standard  value  of  Icm/s  for  the 
gaseous  contaminant 


Table  6-1.  Deposition  rates  and  settling  velocities  for  contamiiumt  species. 

Type 

Deposition  Rate 

Settling  Velocity 

Gaseous 

1.0  cm/s 

0.0  cm/s 

lOpm  particle 

0.0 

0.7 

SOpm  particle 

0.0 

6.3 

lOOpm  particle 

0.0 

52.2 

The  surface  dose  was  computed  on  a  coarse  grid  (20  km  spacing)  covering  most 
of  the  continental  United  States  at  24  and  48  hours  after  the  start  of  the  release.  The 
surface  doses  at  48  hours  for  the  release  of  a  gaseous  material  at  5(X)m  for  each  of  the  3 
days  and  2  times  of  day  are  shown  in  Hguie  6-7  for  the  Illinois  source  and  in  Figure  6-8 
for  the  Four  Comers  source.  The  extent  of  the  surface  pattern,  measured  by  the  area  with 
a  dose  greater  than  10*^  mg/m^  or  the  total  mass  contained  within  that  area  varies  by 
about  a  factor  of  2  excluding  the  2/4/87  Four  Comers  releases.  These  cases  show  a  much 
smaller  plume-like  pattern  but  also  left  the  surface  grid  to  the  west  The  time  of  day  of 
the  release  ((X):00  vs.  12:00)  does  not  significantly  affect  the  amount  of  material  reaching 
the  surface.  Figure  6-9  shows  the  surface  dose  at  48  hours  for  noon  releases  on  2/4/87  at 
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Figure  6-7.  Surface  dose  at  48  hours  after  tee  release  of  a  gaseous  material  at  SOOm  from  tee 
niiuois  source.  Release  start  at  (a)  Ol/07/87:OOkX),  (b)  02A)4/87:00KX).  (c) 
03A)4/87:00KX),  (d)  01A)7/87:12m  (e)  02A)4/87: 12:00,  (0  03/04/87: 12KX). 
Contour  levels  of  10-^,  10^,  10■^  10-X  lO**  and  1  mg/m^. 
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Figure  6-8.  Surface  dose  at  48  hours  after  tee  release  of  a  gaseous  material  at  SOOm  from  the 
Four  Comers  source.  Release  start  at  (a)  01/07/87 :00KX),  (b)  02/04/87:00:00,  (c) 
03/04/87:00:00,  (d)  01/07/87:12:00,  (e)  02/04/87:12:00,  (f)  03/04/87:12:00. 
Contour  levels  of  10-^,  10^,  10*^,  10-^  10**  and  1  mg/m^. 
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Figure  6-9.  Surface  dose  at  48  hours  after  the  release  on  02/04/87  at  12:00  from  the  Illinois 

source.  Release  of  (a)  gas  at  500m,  (b)  gas  at  1000m,  (c)  gas  at  2000m,  (d)  lO^m 
at  500m,  (e)  lO^m  at  lOOOm,  (0 10|tm  at  2000m. 

Contour  levels  of  10-^,  10^,  10-^,  10-2, 10-*  and  1  mg/m^. 
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Figure  6-lO.Surface  dose  at  48  hours  after  the  release  on  02/04/87  at  12:00  from  the  Illinois 
source.  Release  of  (a)  30^m  at  500m,  (b)  3CHim  at  1000m,  (c)  SO^m  at  2000m, 
(d)  lOOvim  at  500m,  (e)  lOO^m  at  1000m,  (f)  lOOtxm  at  200(^. 

Contour  levels  of  10-^,  10^,  10-^  lO**,  10-*  and  1  mg/m^. 
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the  central  lUinoia  site  for  several  release  heights  and  both  a  gaseous  contaminant  and  a 
lOtim  particulate.  Hguie  6-10  shows  the  same  information  for  the  larger  particle 
releases.  The  gaseous  and  10)im  particle  releases  show  nearly  identical  surface  patterns 
and  show  little  variation  with  height  of  release  provided  the  release  is  within  the  mixed 
layer  (<100Qm).  For  releases  above  the  mixed  layer,  the  gaseous  release,  which  is 
dependent  on  turbulent  mixing  to  bring  the  material  to  the  surface,  does  not  impact  the 
surface  while  the  10pm  particles  eventually  Call  into  the  mixed  layer  and  are  mixed  to  the 
surface.  The  largest  particles  (100pm)  are  not  dependent  on  mixing  to  get  to  the  surface 
and  show  little  or  no  effect  of  release  height  or  meteorology. 

The  surface  dose  was  also  computed  on  a  finer  grid  (2  km  spacing)  covering  the 
vicinity  of  the  source  at  4  and  8  hours  after  the  start  of  the  release.  The  surface  dose  at  8 
hours  for  the  release  of  a  gaseous  material  at  20(hn  for  each  of  the  3  days,  and  2  times  of 
day  is  shown  in  Figure  6-11  (Illinois  source)  and  Figure  6-12  (Four  Comers  source).  The 
variations  in  the  extent  of  the  surface  pattern  are  somewhat  less  than  a  factor  of  2  if 
measured  by  the  area  greater  than  10-2  mg/m^  but  greater  than  a  factor  of  2  if  measured  by 
the  mass.  We  should  note  that  the  effect  of  wind  uncertainty  is  evident  in  the  relatively 
wide  spread  of  the  plume.  The  wind  field  is  only  defined  on  a  coarse  grid  and  therefore 
contains  significant  uncertainty  on  the  smaller  scales.  Expected  concentrations  for  a 
particular  event  would  be  higher  than  those  appearing  in  the  figures,  since  these  represent 
an  average  impact  at  each  location  from  all  possible  wind  fields.  The  expected 
concentration  may  display  more  sensitivity  to  source  conditions,  but  we  would  not  expect 
strong  sensitivity  at  the  ranges  and  transport  times  used  in  this  study.  Small-scale 
turbulence  will  effectively  spread  the  plume  at  a  rate  of  order  Itn/s  over  the  first  hour,  so 
variations  in  source  size  of  lO's  of  meters  will  be  negligible,  in  general.  Figures  6-13  and 
6-14  show  the  surface  doses  at  8  hours  for  noon  releases  on  2/4/87  at  the  central  Illinois 
site  for  three  release  heights  and  the  four  release  materials  while  Figures  6-15  and  6-16 
show  the  surface  doses  for  a  night  release.  Close  in,  the  effects  of  release  height  and  time 
of  release  are  more  apparent  for  the  gaseous  and  small  particle  releases.  For  tl^  larger 
particles  release  height  makes  little  difference  but  shows  some  effect  of  the  meteorology. 

In  summary,  the  most  important  source  factors  governing  the  environmental 
prediction  are  the  release  rate  and  duration,  and  the  release  height  The  type  of  mateiial, 
ie.,  gaseous  or  particulate  can  also  be  important  in  determining  the  timescales  of  interest 
for  the  atmospheric  transport  Larga  particles  (diameters  >3()pm)  will  fall  out  relatively 
quickly  near  the  source  and  will  be  less  influenced  by  the  meteorology.  The  mass  release 
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Figure  6-11.  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m  from  the 
Dlinois  source.  Release  start  at  (a)  01/07/87:(n:00,  (b)  02/04/87:00:00,  (c) 
03A)4/87KX):00.  (d)  01/07/87: 12.*00.  (e)  02/04/87: 12.-00,  (f)  03A)4/87: 12:00. 
Contour  levels  of  1(H,  10^^,  lO-^,  lO**,  1  and  10  mg/m^. 


Hgure  6-12.  Surface  dose  at  8  hours  after  the  release  of  a  gaseous  material  at  200m  from  the 
Four  Cornets  source.  Release  start  at  (a)  01/07/87:00:00,  (b)  02/04/87:00:00,  (c) 
03/04/87:00:00,  (d)  01/07/87:12:00,  (e)  02/04/87:12:00,  (f)  03/04/87:12:00. 
Contour  levels  of  10^,  lO-^,  10-2,  jo-i^  i  and  lo  mg/m2. 
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Longitude 

Hguie  6-13.  Surface  dose  at  8  hours  after  die  release  on  02/04/87  at  12:00  from  the  Illinois 

source.  Release  of  (a)  gas  at  20(hn,  (b)  gas  at  1000m,  (c)  gas  at  2000m,  (d)  10^m 
at  200m,  (e)  lO^m  at  1000m,  (f)  IC^  at  2000m. 

Contour  levels  of  10^,  10-^,  lO-^,  10*^  land  10  mg/m^. 


Figure  6-14.  Surface  dose  at  8  hours  after  the  release  on  02/04/87  at  12:00  from  the  Illinois 
source.  Release  of  (a)  30^m  at  200m,  (b)  30^m  at  1000m,  (c)  30^m  at  200(hn, 
(d)  lOOtun  at  200m,  (e)  lOO^m  at  1000m,  (f)  100^  at  200Cto. 

Contour  levels  of  10-*,  lO*^,  10-^,  10-^ ,  1  and  10  mg/m^. 
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gure  6-15.  Surface  dose  at  8  hours  afta  the  release  on  02/04/87  at  00:00  from  the  Illinois 

source.  Release  of  (a)  gas  at  200m,  (b)  gas  at  SOOm,  (c)  gas  at  1000m,  (d)  lOiim 
at  200m,  (e)  lOjim  at  SOOm,  (f)  lO^im  at  1000m. 

Contour  lev^  of  1(H,  lO’^,  lO"*,  lO’^  land  10  mg/m^. 


Longitude 

Figure  6-16.  Surface  dose  at  8  hours  after  the  release  on  02/04/87  at  00:00  from  the  Illinois 
source.  Release  of  (a)  SOpm  at  200m,  (b)  30iim  at  500r4->,  (c)  30iun  at  1000m, 
(d)  lOOiim  at  200m,  (e)  lOOpm  at  SOOm.  (f)  100|im  at  lOOOtn. 

Contour  levels  of  10^,  lO-^,  10-2, 10** ,  1  and  10  mg/m^. 


rate  is  a  scale  factor,  i.e.,  double  the  release  rate  gives  twice  the  concentration 
everywhere,  but  the  quantitative  assessment  of  the  consequences  of  an  actual  release 
clearly  requires  a  reasonably  accurate  definition  of  the  release  mass. 

The  release  height  is  one  of  the  important  parameters  for  the  dispersion 
calculation.  For  long  range  transport  over  time  periods  of  many  hours,  the  key  factor  is 
whether  the  release  is  within  the  planetary  boundary  layer.  A  gaseous  release  above  the 
boundary  layer  will  take  a  long  time  to  diffu%  vertically  and  may  not  reach  the  surface 
unless  it  encounters  deep  mixing,  e.g.,  in  the  vicinity  of  convective  clouds.  Conversely,  a 
release  within  the  turbulent  boundary  layer  will  be  quickly  mixed  throughout  the  depth  of 
the  layer  and  is  therefore  only  sensitive  to  the  initial  height  very  close  to  the  source.  This 
region  of  sensitivity  extends  a  few  release  heights  downstream,  since  the  vertical 
turbulent  velocity  fluctuations  are  generally  some  fraction  of  the  mean  wind  speed,  and  is 
therefore  of  the  order  of  a  few  hundred  meters  to  a  few  kilometers.  Beyond  this  range, 
the  contaminant  will  be  well-mixed  under  normal  daytime  conditions. 

For  nocturnal  releases,  the  short  term  prediction  is  critically  dependent  on  the 
description  of  the  boundary  layer.  Details  of  the  boundary  layer  are  not  generally 
available,  so  our  knowledge  of  the  vertical  structure  is  rarely  more  reliable  than  a  few 
hundred  meters.  It  seems  inappropriate  for  a  source  description  to  be  more  accurate  than 
this. 
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SECTION  7 


CONCLUDING  REMARKS 

The  focus  of  this  report  has  been  the  development  and  improvement  of  techniques 
for  predicting  the  late-time  dispersion  of  nuclear  clouds  in  the  atmosphere,  with  particular 
emphasis  on  the  characterization  of  the  uncertainties  associated  with  such  a  prediction. 
Our  knowledge  of  the  state  of  the  atmosphere,  either  in  real  time  or  at  some  future  time, 
is  inevitably  limited,  and  this  inability  to  specify  the  complete  wind  field  generally  leads 
to  large  uncertainty  in  any  forecast  of  dispersion.  In  many  instances,  the  range  of  the 
prediction  errors  is  comparable  to  or  larger  than  the  predicted  value  so  it  is  very  important 
to  be  able  to  describe  the  statistical  distribution  of  these  errors.  The  errors  cannot  be 
reduced  without  more  detailed  information  about  the  atmosphere,  and  since  this  is  a 
formidable  or  even  impossible  (in  the  case  of  future  scenarios)  taric,  it  is  important  that  a 
probabilistic  prediction  cs^ability  be  developed. 

The  basis  of  the  probabilistic  prediction  methodology  reported  herein  is  the 
second-order  turbulence  closure  framework  as  mbodied  in  the  SCIPUFF  model,  which 
uses  a  prognostic  equation  for  the  variance  in  the  predicted  concentration  level  to  assess 
the  level  of  uncertainty  in  any  predicted  quantity.  The  model  has  been  extended  to 
describe  large  scale  atmospheric  dispersion,  using  the  GASP  velocity  fluctuation  data  of 
Nastrom  and  Gage  (1985)  and  the  velocity  spectra  distribution  concepts  of  Gifford 
(1988).  Although  limited  in  their  representation  of  atmospheric  variability,  since  the 
GASP  measurements  were  taken  by  commercial  aircraft  on  international  flights,  these 
data  provide  an  irutial  foundation  on  which  to  build  a  rational  description  of  long  range 
di^)ersiorL  Furthermore,  Bauer  (1983)  concludes  that  the  dispersion  characteristics  of  the 
atmosphere  do  not  vary  marioedly  with  altitude,  so  that  the  aircraft  cruise  altitude  data  can 
be  used  as  a  general  basis  for  horizontal  dispersiorL 

In  addition  to  the  improvement  in  the  turbulence  specification,  the  Gaussian  puff 
description  in  SCIPUFF  has  been  extended  to  include  a  completely  general  spatial 
second-moment  evolution  equation.  This  allows  individual  puffs  to  distort  under  the 
action  of  wind  shear  without  artificially  diffusing  the  concentration  field  or  allowing 
neighboring  puffs  to  separate  from  one  another,  lea^g  gaps  in  the  initially  continuous 
cloud.  The  crq)ability  of  the  new  scheme  to  describe  the  persistent  action  of  shear  on  an 


isolated  cloud  was  demonstrated  in  the  'Dust-OfT  model  comparison  exercise,  where 
contours  of  die  concentration  field  after  8  hours  transport  showed  a  smoothly  continuous 
distribution  extending  over  lOOCKon  in  tte  horizontal.  The  model  comparison  statistics 
were  not  compiled  as  part  of  this  study,  but  the  treatment  of  wind  shear  was  one  of  the 
most  sensitive  aspects  of  the  various  model  approaches.  Maintaining  a  strongly  sheared 
cloud  without  artificial  diffusion  is  a  difficult  numerical  requirement,  but  is  a  prevalent 
feature  of  high  altitude  dispersion. 

Model  evaluation  is  always  a  vital  part  of  the  development  of  a  predictive 
capability.  Held  data  on  atmospheric  dispersion  over  scale  of  hundreds  to  thousands  of 
kilometers  is  relatively  sparse,  however,  making  it  difficult  to  quantitatively  assess  the 
accuracy  of  a  model.  The  Across  North  America  Tracer  Experiment  (ANATEX) 
provides  a  uniquely  extensive  data  base  for  model  evaluation  at  continental  scale,  since 
data  were  taken  on  an  extensive  network  over  a  period  of  three  months.  This  experiment, 
sponsored  by  NOAA,  used  specialized  tracers  detectable  at  very  low  concentrations  with 
periodic  surface  releases  throughout  the  experiment.  In  addition,  NOAA  provided 
detailed  wind  fields  for  model  evaluation  purposes.  SCIPUFF  was  used  to  simulate  the 
entire  experimental  period  and  the  statistical  distributions  of  concentrations  were 
compared  with  the  observations.  The  model  showed  good  agreement  with  the  observed 
24-hour  average  concentrations  as  a  function  of  distance  from  the  source,  and  also 
showed  good  agreement  with  the  short-term  aircraft  samples  obtained  on  flights  through 
the  plume.  SCIPUFF  also  showed  some  skill  in  predicting  the  observed  distributions 
using  only  a  single  statistical  average  description  of  the  wind  field  over  the  entire 
continent  for  the  three  month  period.  The  evaluation  exercise  demonstrated  the  capability 
of  the  model  in  predicting  concentration  distributions  for  both  long-term  and  short-term 
averages,  and  also  indicates  some  skill  using  very  limited  wind  data. 

Several  studies  are  reported  in  the  preceding  sections  that  demonstrate  the 
feasibility  of  a  quantitative  treatment  of  uncertainty  in  atmospheric  dispersion 
predictions.  The  ANATEX  evaluation  exercise  provides  some  real  data  for  comparison 
with  SCIPUFF  concentration  predictions,  but  the  probabilistic  framework  clearly  has 
wider  applicability  and  requires  further  validation.  The  general  basis  for  the  SCIPUFF 
prediction  of  turbulent  boundary  layer  dispersion  now  seems  to  be  established  on  scales 
ranging  from  several  kilometers  to  thousands  of  kilometers,  but  there  are  a  number  of 
areas  where  further  development  is  required.  For  example,  our  understanding  of  the 
variation  in  the  one-point  probability  distribution  function  with  averaging  time  and 
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transport  distance  could  be  improved.  The  ANATEX  results  suggest  that  the  clipped- 
normal  distribution  is  valid  for  short  term  samples,  but  long  term  averages  are  better 
approximated  as  log-normal.  A  second  area  of  research  is  the  high  altitude  shear- 
dominated  regime,  where  the  vertical  diffusion  is  suppressed  by  the  atmospheric  stability. 
There  is  little  data  available  from  this  regime,  and  our  extension  of  the  prediction 
methodology  from  the  relatively  well-understood  turbulent  lower  atmosphere  is  not  based 
on  solid  understanding. 

In  summary,  therefore,  we  have  demonstrated  a  practical  methodology  for  the 
prediction  of  late-time  dispersion  in  a  probabilistic  framework.  The  SCIPUFF  model  is 
capable  of  providing  routine  forecasts  of  concentrations  from  arbitrary  sources  together 
with  a  statistical  assessment  of  the  possible  impact  The  model  has  been  validated  on  a 
wide  range  of  horizontal  scales,  using  surface  tracer  data,  and  shows  good  agreement 
with  concentration  measurements.  The  probabilistic  framework  should  be  extended  and 
improved  for  use  in  general  dispersion  calculations,  and  the  evaluation  studies  should  be 
continued. 
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ABSTRACT 

A  lonB>ranfe  transport  model  based  on  turbulence  closure  concepts  is  described.  The  model  extends  the 
description  of  planetary  boundary  layer  turbulent  diffusion  to  the  laiier  scales  and  uses  statistical  wind  information 
to  predict  conuminant  dispersion.  The  model  also  contains  a  prediction  of  the  sutistical  fluauations  in  the 
tracer  concentration  resulting  from  the  unresolved  velocity  fluauations.  The  dispersion  calculation  is  made  by 
means  of  a  Lagrangian  puff  represenution.  allowing  the  use  of  time-dependent  three-dimensional  flow  fields. 
Predictions  of  the  ANATEX  ( Across  North  America  Tracer  Experiment )  releases  are  compared  with  observauons. 

Both  24-h  average  surftoe  and  short-term  aircraft  sampler  concentrauons  are  calculated  using  the  h«li-iesoluuon 
wind  fields  from  the  NMC  Nested  Grid  Model.  The  sutistical  prcdicuon  is  also  tested  using  long-term  average 
wind  dau. 

Sutistical  uncertainty  in  the  predictions,  due  to  the  unresolved  wind  fluauauons.  is  found  to  be  small  for 
the  24-h  average  surface  concentrations  obtained  with  the  h^resolution  winds  but  is  very  significant  for  the 
short-term  aircraft  sampler  concentrations.  A  clipped  normal  probability  distribuuon  provides  a  reasonably 
good  descnpuon  of  the  overall  cumulative  distribution  of  the  aircraft  sampler  concentrations.  A  reasonably 
good  description  of  the  24-h  surface  concentrations  is  also  obtained  using  only  the  long-term  average  wind 
sutisucs  and  a  lognormal  probability  distribuuon  for  the  concentiauon  values. 


1.  iBtrodactioa 

The  atmospheric  transport  of  pollutants  over  con¬ 
tinental  scales  is  an  important  environmental  process 
and  has  received  considerable  attention  from  modelers 
and  experimentalists.  With  national  policy  on  issues 
such  as  acid  rain  being  influenced  by  the  prediaions 
of  long-range  transport  models,  there  is  a  strong  de¬ 
mand  for  verification  of  these  models.  As  with  other 
dispersion  model  evaluation  studies  at  smaller  scales, 
the  comparison  between  field  observations  and  model 
prediaions  is  far  from  straightforward.  The  inability 
to  determine  the  state  of  the  atmosphere  leaves  a  large 
component  of  uncertainty  in  the  prediction.  It  is  gen¬ 
erally  very  difficult  to  distinguish  model  error  from 
maeorological  input  uncertainty,  and  consequently,  it 
is  impossible  to  reliably  assess  model  accuracy.  Various 
statistical  measures  have  been  recommended  for  model 
evaluation  (Fox  1981).  but  there  are  no  standard 
methods  designed  to  address  the  uncertainty  in  a 
quantitative  way.  As  Venkatram  ( 1982)  points  out. 
the  statistical  ensemble  represented  by  any  model  is 
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defined  in  terms  of  the  model  input  data,  which  makes 
comparison  between  different  models  difficult  if  they 
use  different  input. 

Recognizing  the  stochastic  nature  of  atmospheric 
dispersion,  we  have  previously  proposed  that  a  quan¬ 
titative  prediaion  of  the  uncertainty  be  included,  al¬ 
lowing  statistical  significance  to  be  assessed  ( Lewellen 
and  Sykes  1989).  The  uncenainty  is  inherent  in  the 
sense  that  the  atmosphere  can  never  be  complaely  de¬ 
scribed.  and  effort  should  be  expended  to  d^  direaly 
with  the  random  component  rather  than  toward  elim¬ 
inating  it.  We  have  therefore  developed  models  for  at¬ 
mospheric  dispersion,  based  on  second-order  turbu¬ 
lence  closure,  that  include  a  prediction  of  the  concen¬ 
tration  fluauation  variance  and  thus  provide  a  means 
of  testing  consistency  between  observations  and  pre¬ 
dictions.  These  models  were  designed  to  predia  plume 
dispersion  out  to  ranges  of  about  SO  km  and  were  ex- 
tenavely  compared  with  the  Elearic  Powa  Research 
institute  ( EPRl )  plume  model  validation  experiments 
( Lewellen  et  al.  1988). 

The  major  advantage  of  the  use  of  turbulence  closure 
in  the  dispersion  prediaion  is  the  ability  to  relate  the 
evolution  of  the  scalar  concentration  field  to  measur¬ 
able  velocity  statistics.  This  approach  has  been  initiated 
by  Gifford  ( 1988)  for  long-range  transport  prediaions 
using  a  random  walk  model  and  is  presented  here  using 
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the  second-order  closure  framework.  The  dispersion 
calcuiauon  is  extended  to  include  a  prediaion  of  the 
statistical  fluauation  about  the  mean  value.  This  pro¬ 
vides  a  method  of  ^xounting  for  the  effects  of  averaging 
time  on  the  observed  concentration  as  a  function  of 
the  turbulent  wind  input.  In  this  context,  turbulence 
denotes  the  unresolved  component  of  the  wind  field 
and  therefore  depends  strongly  on  the  wind  data  used 
in  the  model  prediction.  The  statistical  framework 
forces  explicit  consideration  of  the  effects  of  averaging 
time  and  eddy  scales,  in  contrast  to  the  implicit  as¬ 
sumptions  inherent  in  the  choice  of  diffusion  param¬ 
eters. 

In  this  paper,  we  describe  the  adaptation  of  our 
Gaussian  puff  model  to  the  long-range  dispersion 
problem  and  compare  the  model  results  with  the  tracer 
observations  from  the  Across  North  America  Tracer 
Experiment  (ANATEX)  study  (Draxler  et  al.  1991 ). 
ANATEX  was  conducted  in  1987  under  the  auspices 
of  the  National  Oceanic  and  Atmospheric  Adminis¬ 
tration  ( NOAA )  and  provides  a  high  quality,  extensive 
database  for  dispersion  over  thousands  of  kilometers. 
The  daubase  is  large  enough  for  statistical  comparisons 
to  be  significant  and  appears  to  be  the  most  complete 
study  to  date. 

The  next  section  briefly  describes  the  ANATEX  data 
employed  in  our  model  evaluation  study.  Section  3 
details  the  puff  model  as  extended  to  the  continental- 
scale  problem  and  discusses  the  model  comparison 
statistics.  Section  4  contains  the  results.  Concluding 
remarks  are  conuined  in  seaion  S. 

2.  ANATEX  data 

a.  ANATEX  experiment 

The  ANATEX  experiment  is  fully  described  in  a 
senes  of  repons  ( Draxler  and  Heffter  1989;  Stundler 
and  Draxler  1989;  Heffter  and  Draxler  1989)  and 
summanzed  by  Draxler  et  al.  ( 1 99 1 ) .  so  only  the  most 
important  features  will  be  outlined  here.  The  experi¬ 
ment  was  designed  to  provide  a  database  for  long-range 
tianspon  models  and  used  an  inen  tracer  detectable 
at  very  low  concentrations  to  track  the  dispersion  across 
the  continental  United  States.  Small  releases  of  per- 
fluorocarbon  tracers  were  made  from  Glasgow.  Mon¬ 
tana.  and  from  St.  Cloud.  Minnesou.  and  measured 
on  a  network  of  ground  samplers  extending  eastward 
from  the  sources  to  the  Atlantic  coast  and  southward 
to  the  Gulf  of  Mexico.  The  release  locauons  and  surface 
sampler  network  are  shown  in  Fig.  I . 

Three  tracers  were  employed,  with  detection  limits 
below  1  fl  r'  =  10"‘*  1 1*',  The  experiment  lasted  3 
months,  from  5  January  through  30  March  1987.  with 
releases  made  for  3  h  every  2.5  days  from  the  two  source 
locations.  Alternate  releases  were  made  in  daytime  and 
nocturnal  conditions.  Perfluoro-orthodimethylcyclo- 


Fic.  I.  ANATEX  pnnwry  ground  level  sampling  network. 


hexane  (oPDCH)  and  Perfluoromethylcyclohexanc 
( PMCH )  were  released  from  St.  Cloud,  with  PMCH 
being  released  every  5  days  only  under  daytime  con¬ 
ditions  to  provide  some  measure  of  independence. 
Perfluorotrimethylcyclohexane  ( PTCH )  was  released 
only  from  Glasgow.  There  was  some  small  conumi- 
nation  of  the  oPDCH  release  by  PMCH,  but  this  seems 
to  be  adequately  accounted  for  in  the  release  dau.  The 
tracer  technology  is  extensively  described  by  Draxler 
and  Heffter  ( 1989).  and  the  reader  is  referred  to  the 
ANATEX  reports  for  further  details. 

b.  Ground  samplers 

The  ground  sampler  network,  illustrated  in  Fig.  1. 
shows  a  toul  of  77  collection  sites.  Each  sampler  col¬ 
lected  data  for  24-h  periods,  which  were  then  analyzed 
in  the  laboratory  to  obtain  the  average  concentrations 
of  all  three  tracers  over  the  penod.  Samples  were  an¬ 
alyzed  in  two  independent  laboratones  to  provide  es¬ 
timates  of  reliability,  and  several  collocated  samples 
were  also  uken  to  determine  repeaubility.  Estimates 
of  the  background  concentrations  of  the  three  tracers 
were  made  from  the  distribution  of  analyzed  concen¬ 
trations.  and  all  sampler  concentrations  are  reported 
as  excess  over  the  background  levels.  There  is  some 
uncertainty  in  the  estimate,  and  excess  values  close  to 
the  background  value  are  not  reliable.  The  background 
values  for  all  three  tracers  are  shown  in  Table  I . 

In  general,  the  sampler  data  analysis  provided  reli¬ 
able  concentrations  from  a  very  large  fraction  of  the 
total  data  colleaed.  Roughly  84%  of  the  total  possible 
number  of  measurements  were  archived  as  reliable 
dau.  giving  5400  concentrations  for  each  of  the  three 
tracers.  This  is  a  remarkably  high  return,  and  the 
quantity  and  quality  of  the  dau  make  significant  su- 
tistical  comparison  possible  using  this  dataset. 
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Tasle  I.  Ambiciit  tracer  background  concentrations. 
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c.  Aircraft  data 

In  addition  to  the  surface  sampler  network,  three 
aircraft  made  direct  sampling  flights  through  the 
plumes  at  ranges  out  to  about  450  km  from  the  source. 
The  aircraft  attempted  to  fly  transverse  paths  through 
the  plume  at  different  altitudes,  using  observed  and 
forecast  wind  fields  to  estimate  the  plume  location.  The 
flights  were  reasonably  successful,  with  roughly  40%  of 
the  observations  producing  nonzero  concentrations, 
indicating  the  presence  of  a  plume.  Flights  were  made 
at  altitudes  up  to  about  2500  m.  with  data  being  sam¬ 
pled  over  6-min  periods  ( 1 2  min  in  some  cases).  Using 
a  typical  airspeed  of  about  55  m  s  ~ .  the  effective  path- 
length  was  about  20  km  ( or  40  km). 

In  total,  the  aircraft  provided  about  1 300  sampler 
measurements  from  a  total  of  50  flights  from  Glasgow, 
and  22  flights  from  St.  Goud.  The  short  time-average 
path  integrals  are  an  interesting  complement  to  the 
fixed  24-h  averages  ftom  the  ground  samplers.  Unfor¬ 
tunately,  the  limited  range  of  the  aircraft  ^ta  prevents 
much  direct  comparison  between  the  two.  The  aircraft 
dau  does  present  a  serious  challenge  for  model  pre¬ 
diction.  however,  and  the  predictability  of  these  near- 
instantaneous  samples  will  be  examined  below. 

3.  Model  descriptkw 

a.  Puff  dispersion  model 

SCIPUFF  f  for  second-order  closure  integrated  puff) 
was  originally  developed  as  a  shon-range  dispersion 
model  and  tested  against  the  EPRl  plume  model  de¬ 
velopment  and  validation  experiments  at  the  Kincaid. 
Illinois,  and  Bull  Run.  Tennessee,  power  plants  ( Sykes 
et  al.  1988).  The  model  was  derived  from  simplifica¬ 
tions  of  the  full  second-order  turbulence  closure  equa¬ 
tions  using  a  Gaussian  puff  ftamework.  The  Lagrangian 
puff  methodology  allows  a  general  computation  of 
time-dependent  effects  and  spatial  inhomogeneity.  The 
model  has  been  modified  considerably  to  allow  appli¬ 
cation  to  long-range  transport  problems,  and  the  full 
equations  for  a  passive  tracer  are  presented  here.  SCI¬ 
PUFF  represents  the  scalar  concentration  distribution 
as  a  colleaion  of  overlapping  Gaussian  puffs:  that  is, 

C(x,t)=  2  Caix.l).  (1) 

o 


where,  using  sundard  Cartesian  tensor  notation  and 
the  summation  convention. 
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1/211 


1/2 


X  exp 
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Here  Q„  is  the  mass  assoaated  with  puff  a.  .v,„  is  the 
puff  centroid.  a„o  is  the  second-moment  tensor,  and 
lia.il  is  the  determinant  of  a,„.  The  formal  definition 
of  the  moments  is 


0.  = 

[cjv. 

(3) 

0.9^ ia  ~~ 

J  CaX,dr. 

(4) 

II 

CjX,  -  X.a)(Xj  -  Xja)dl'. 

Jv 

(5) 

This  is  a  generalization  of  the  sundard  puff  represen- 
ution.  which  assumes  that  a  is  diagonal  (e.g..  Bass 
1980:  Sykes  etal.  1988). 

The  model  is  advanced  in  time  by  solving  ordinary 
differential  «)uations  for  the  puff  moments  ( 3 )-( 5 ). 
The  derivation  of  these  equations  follows  the  plume 
aiuiysis  of  Sykes  et  al.  ( 1986)  and  will  be  described 
only  briefly. 

Conservation  of  scalar  mass  implies 


(6) 


Consideration  of  the  moment  equations  yields 


M,(X..  1) 


(7) 


and 


da,,, 

dt 


du,  du, 

-  F.JO  +  ^  ^  • 


(8) 


where  u,  is  the  ambient  mean  wind  velocity,  and 


0-f.A.  =  j  l<X,  -  X,a)U',C'a  +  (X;  -  Xja)U',Ca]dV. 

(9) 

Here  u',c'a  represents  the  turbulent  flux  of  the  scalar 
associated  with  puff  a:  sundard  Reynolds-averaging 
noution  is  used,  with  the  overbar  representing  the 
mean  value  and  prime  denoting  fluauation  from  the 
mean.  Thus.  ( 7 )  represents  translation  of  the  centroid 
by  the  mean  wind,  while  ( 8 )  describes  spreading  of  the 
puff  by  turbulent  diffusion  and  distortion  by  wind 
shear. 

The  evolution  of  the  turbulent  diffusion  terms  is  ob¬ 
tained  from  the  second-order  closure  result  of  Sykes  et 
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al.  ( 1986)  but  uses  a  simplified  treatment  of  the  large* 
scale  turtnilence.  yielding 


and 


where 


at  A« 

^  G.  -  /I  fjj. 

ut  Jq  A^- 


dGa  -I-  <7« 

— G.. 

ut  Om  a\|> 
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(10) 

(11) 
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Here  u'* ,  r'* ,  m-'*  are  the  ambient  velocity  variances, 
gy  *  {/'*  +  v'-.  and  Ah  is  the  horizontal  length  scale 
of  the  velocity  fluctuations.  The  parameter  Ar  is  a  ver* 
tical  flunuation  scale,  and  g]  «  3m’'*  is  a  measure  of 
the  energy  associated  with  the  vertical  motions.  The 
gravitational  acceleration  is  g.  To  is  a  reference  tem¬ 
perature  for  the  Boussinesq  approximation,  and  0  is 
the  ambient  potential  temperature.  The  ambient  ver¬ 
tical  heat  flux  is  w'd',  and  G«  represents  the  effects  of 
buoyancy  on  the  vertical  difiusivity.  The  empirical  tur¬ 
bulence  closure  constant.  A,  takes  the  i^ue  0.7S 
(Lewellen  1977).  We  assume  that  »  0  for  /  #  j, 
although  the  effects  of  velocity  cross  correlations  could 
be  included  if  required.  We  do  not  expect  significant 
eftcts  from  such  correlations  in  dispersion  on  the  scales 
of  interest  here.  The  specification  of  the  ambient  wind, 
temperature,  and  turbulence  information  will  be  pre¬ 
sented  in  the  next  section. 

To  account  for  the  inversion-capped  mixing,  the 
vertical  diffusivity  is  set  equal  to  zero  for  puffs 
within,  the  mixed  layer  with  C;  >  0.8r, ,  where 
*  vTii  and  z,  is  the  mixed-layer  depth.  This  simply 
turns  off  any  fiirther  vertical  diffusion  when  the  puff  is 
mixed  across  the  entire  boundary  layer. 

A  valuaUe  feature  of  the  secon^rder  closure  model 
for  the  turbulent  diffusion  is  the  quantitative  link  that 
it  provides  with  measurable  turbulence  correlations. 
The  velocity  variances  are  not  necessarily  easily  mea¬ 
sured.  but  they  can  be  estimated  at  various  length  scales 
to  give  a  direct  prediction  of  scalar  diffusion. 

Having  presented  the  model  equations  for  scalar 
transport  and  diffusion,  it  remains  to  define  the  scalar 
source.  Details  of  the  source  specifications  are  given  in 
appendix  A.  but  a  brief  overview  is  given  here.  Since 
the  dispersion  calculation  time  step  is  15  min  and  the 
source  diameter  is  small  ( less  than  a  meter),  the  spatial 
evolution  of  the  source  plume  must  be  consider^  be¬ 
fore  initializing  an  appropriate  puff.  Essentially,  a 


simplified  one-dimensional  (or  plume)  version  of 
SCI-PUFF  is  used  to  calculate  mass,  centroid,  and 
second  moments  using  a  number  of  small  time  steps 
for  a  total  time  r,.  This  scheme  also  allows  us  to 
compute  the  rapid  evolution  of  the  fluctuation 
variance  c'*. 

Equation  ( 8 )  shows  that  puff  size,  shape,  and  ori¬ 
entation  change  through  the  aaion  of  turbulent  dif¬ 
fusion  and  sheanng  motion  due  to  nonuniformities  in 
the  wind  field.  These  tend  to  increase  the  puff  size; 
diffusion  increases  the  diagonal  moments  and  hence 
increases  the  puff  volume,  while  the  velocity  shear 
terms  preserve  volume  but  cause  elongation  in  the  di¬ 
rection  of  the  wind  shear  (and  contraction  normal  to 
that).  At  some  point,  then,  the  puff  size  will  grow  too 
large  for  the  local  specification  of  the  wind  field  to  be 
accurate  over  the  entire  puff.  This  was  not  a  factor  in 
the  limited-range  modeling  of  Sykes  et  al.  ( 1988)  but 
will  be  important  in  long-range  transport.  Therefore, 
a  puff-splitting  algorithm  based  on  the  diagnnai  com¬ 
ponents  of  the  moment  tensor  a  has  been  imple¬ 
mented. 

When  the  moment  in  a  particular  coordinate  direc¬ 
tion  exceeds  the  critical  size  specified  for  that  direction, 
the  puff  is  split  into  two  smaller,  overlapping  puffs. 
The  critical  size  is  related  to  the  scale  of  the  specified 
mean  velocity  variation  in  each  direction,  typic^y  the 
grid  length  for  wind  fields  from  a  numerical  model 
output.  The  splitting  algorithm  conserves  ail  the  mo¬ 
ments  of  the  distribution,  eliminating  artificial  diffu¬ 
sion.  and  there  is  sufficient  overlap  of  the  newly  created 
pair  of  puffs  to  prevent  significant  change  in  local  con” 
centrations. 

Puib  are  split  in  the  vertical  direaion  as  well  as  in 
the  horizontal,  and  special  care  must  be  used  near  the 
surface  and  near  the  inversion.  The  algorithm  used  in 
SCIPUFF  reflects  puffs  at  the  surface  and  also  at  the 
inversion  for  any  puff  splitting  from  within  the  mixed 
layer.  This  essentially  treats  the  inversion  as  a  rigid 
boundary,  except  that  the  inversion  moves  vertically 
and  the  reflection  condition  is  applied  only  locally  to 
determine  the  location  of  the  newly  created  split  puff. 
A  puff  will  cross  the  inversion  if  the  inversion  falls 
below  its  current  location,  leaving  it  in  the  low-tur¬ 
bulence  region  above  the  mixed  layer. 

In  addition  to  the  splitting  algorithm.  SCIPUFF  also 
incorporates  a  merging  routine  for  puffs  that  sufficiently 
overlap.  Merging  is  required  to  track  diffusion  from  a 
small  source  where  a  plume  of  small  puffs  is  released. 
In  this  case,  the  puffs  will  eventually  grow  and  overlap 
one  another  so  that  the  plume  can  be  adequately  rep¬ 
resented  by  a  smaller  number  of  puffs.  The  merging 
algorithm  maintains  an  efficient  calculation  of  the 
subgrid  concentration  field  and  also  prevents  the  split¬ 
ting  algorithm  from  generating  puffs  in  exponentially 
growing  numbers.  When  new  puffs  are  created  in  re¬ 
gions  that  already  contain  similar  puffs,  they  will  be 
merged  with  the  existing  puffs.  The  merging  procedure. 
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like  the  splitting  algorithm,  maintains  all  the  moments 
of  the  puff  distributions. 

b.  Concentration  fluctuation  variance 

One  of  the  most  imponant  aspects  of  our  modeling 
of  scalar  dispersion  has  been  the  derivation  of  model 
equations  for  the  concentration  fluauation  variance 
7^.  The  turbulent  dispersion  of  a  scalar  is  inherently 
random  and  a  quantitative  measure  of  the  uncertainty 
in  any  observable  concentration  is  almost  essential  for 
any  meaningful  interpretation  of  model  prediaions  of 
atmospheric  observations.  The  plume  models  described 
by  Sykes  et  al.  ( 1984,  1986)  contain  a  variance  pre¬ 
diction  that  has  been  adapted  to  the  context  of  long- 
range  transit. 

The  variance  prediction  is  obtained  from  the  con¬ 
servation  equation  for  mean-square  concentration.  It 
is  easily  shown  that  c~  is  adveaed.  diffused,  and  dis¬ 
sipated  only  by  molecular  diffusivity.  We  a^ume  that 
the  advection  and  turbulent  diffusion  of  c~  is  similar 
to  that  of  the  mean  concentration,  and  represent  the 
c~  field  as  a  sum  of  Gaussians  with  the  same  shape  as 
the  mean.  Thus,  an  additional  puff  quantity.  5.,  is  as¬ 
sociated  with  each  puff,  where 

(15) 

Ji; 

The  evolution  of  is  given  by 


The  representation  of  the  overlap  integral  in  (18) 
relies  on  the  assumption  that  the  integrated  c '  asso¬ 
ciated  with  puff  a  can  be  wntten 


Thi.^  is  a  natural  partition  of  the  total  c*  between  the 
puffs  but  is  somewhat  arbitrary,  and  perhaps  other  par¬ 
titions  could  be  imagined. 

Having  defined  I'.,  all  that  remains  is  the  specifi¬ 
cation  of  the  dissipation  time  scale,  r^.  As  in  our  earlier 
modeling,  is  obtained  from  a  prediaed  length  scale 
of  the  fluctuations.  In  the  long-range  transpon  prob¬ 
lems.  however,  the  difference  in  charaaer  between 
horizontal  and  vertical  motions  must  be  dealt  with. 

The  horizontal  fluctuation  scale  is  modeled  in  a 
similar  way  to  the  plume  analysis  of  Sykes  et  al.  ( 1986 ). 
Thus. 

d  0.25q, .  ///  <  All 

(20) 

ot  [O.lfl,.  /h>A„. 


where 

q;  =(«'*-!-  w'^)  min 
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where  c.  represents  the  volume-integrated  dissipation 
associated  with  puff  a. 

Models  for  fluctuation  dissipation  are  generally 
written  in  the  form,  c'-lr^,  that  is.  the  fluctuation  vari¬ 
ance  divided  by  a  dissipation  time  scale.  The  first  task 
in  the  specification  of  is  therefore  a  derivation  of  the 
variance  from  S„.  By  definition,  c'*  -  c'  -  c'.  The 
variance  can  be  obtained  from  this  relation,  but  the 
calculation  of c'  involves  an  interaaion  between  over¬ 
lapping  puffs.  Since  c=Z„  c^,  we  have 

C- =  "Z  C„c,i.  (17) 

a.d 

We  define  c'jdV,  as  the  integrated  variance 

associated  with  puff  a.  and  approximate 

i;  =  -  V  (18) 

a 

where  is  the  overlap  integral  between  puff  a  and 
puff  d:  that  is. 

/..d  =  J««(x)gfl(x)dl’.  (19) 

andg,  =  Ca/Qa  from  (2). 


This  represents  an  inertial  range  scaling  of  the  ve¬ 
locity  variances,  with  the  vertical  component  included 
to  account  for  the  increased  energy  in  conveaive  con¬ 
ditions. 

We  then  model 

T:i'i  =  2bsj^.  (22) 

in 

where  6  =  0.125  and  5  =  1.8  are  closure  model  con¬ 
stants  ( Lewellen  1977). 

The  total  dissipation  must  also  account  for  vertical 
motions  since  the  diffusion  through  the  depth  of  the 
convective  layer  will  involve  significant  dissipation  of 
concentration  fluctuations.  In  the  context  of  long-range 
transport,  it  is  assumed  that  this  is  a  rapid  process:  that 
is.  the  pollutant  quickly  becomes  well-mixed  through 
the  boundary  layer  and  does  not  produce  fluctuations. 
Generally,  the  time  scales  associated  with  vertical  mo¬ 
tions  are  expeaed  to  be  rapid,  although  they  will  not 
remove  the  large-scale  horizontal  fluauations.  There¬ 
fore.  a  dissipation  term  has  been  included  to  reduce 
the  total  S„  by  the  appropriate  amount  to  match  the 
reduaion  of  the  mean  concentration  by  vertical  dif¬ 
fusion.  The  appropnate  dissipation  rates  are 
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and 


(24) 


The  full  model  for  5.  is  therefore 

dt  "  T,  ’ 

where  '  =  t,//  +  .  The  vertical  puff  scale  is  rep¬ 

resented  by  iv,  which  is  not  the  same  as  a..  The  indi¬ 
vidual  puff  spread  depends  on  the  splitting  algorithm, 
which  will  reduce  a.  to  maintain  numerical  resolution. 


the  layer.  Under  our  assumption  of  rapid  mixing  for 
long-range  transport,  the  detailed  evolution  of  the 
mixed-layer  profile  is  not  important.  Difliculties  are 
avoided  with  the  splitting  procedure  by  specifying  con¬ 
stant  values  that  are  representative  of  the  entire  layer. 
Thus, 


0.1  ui  +  0.4^^,  :  < 

lO'^m^s"^,  :>z,, 

max(2ui  +  0.4wi.  u't),  z  k  z, 
w'r.  z>z,. 


c.  Meteorology 

The  meteorological  inputs  necessary  for  the  disper¬ 
sion  calculation  are  the  mean  wind  field  and  the  tur¬ 
bulent  velocity  correlations.  SCIPUFF  does  not  contain 
any  meteorological  data-interpolation  modules,  so  the 
meteorology  must  be  specific.  Fortunately,  the  Na¬ 
tional  Meteorological  Center's  Nested  Grid  Model 
( NGM )  output  was  archived  every  2  h  throughout  the 
ANATEX  experiment  and  is  available  on  tape  from 
the  NOAA  Air  Resources  Laboratory.  This  high-res¬ 
olution  weather  prediction  model  was  run  in  12-h 
forecast  mode  for  the  period,  giving  short-range  fore¬ 
casts  on  a  horizontal  grid  of  roughly  90  km.  For  use 
in  SCIPUFF,  the  fields  were  interpolated  onto  a  three- 
dimensional  grid  with  1**  horizontal  (latitude,  longi¬ 
tude)  resolution  and  nonuniform  vertical  resolution 
as  indicated  in  Table  2. 

The  boundary-layer  turbulence  was  deduced  from 
the  NGM  predictions  of  surface  stress,  ui,  where 
is  the  surface  friction  velocitv,  heat  flux.  Ho,  and  mixing 
layer  depth,  r, .  as  descri'  below.  Single  mixed-layer 
values  of  the  turbule.  ..  isity  and  length  scale  for  a 

given  latitude  and  lo.<K>i JO  tvere  used  because  the  puff¬ 
splitting  algorithm  tends  to  produce  a  nonuniform 
concentration  profile  if  the  mixed-layer  contains  several 
vertical  grid  lengths  and  the  diffusivity  varies  across 


Table  2.  Vertical  grid  used  for  SCIPUFF  meteorology  (m). 


Uvei 

Height 

1 

0.0 

150.0 

3 

300.0 

4 

485.7 

5 

715.5 

6 

999.6 

7 

1350.4 

8 

1783.2 

9 

2316.9 

10 

2974.2 

II 

3783.2 

i: 

4777.9 

13 

6000.0 

14 

7500.0 

X, 


'0.25r,.  r<r, 
10  m  z>z,. 


where  w*  =  [{f/  To)HoZ,  ]  is  the  convective  velocity 

scale  and  the  free  atmosphere  horizontal  velocity  vari¬ 
ance  is  discussed  in  the  next  section.  Here  Ho  is 
assumed  positive  in  the  definition  of  w, .  If  //o  is  neg¬ 
ative.  w,  is  set  to  zero.  The  free  atmosphere  values  for 
the  vertical  velocity  variance  and  scale  are  somewhat 
arbitrary,  although  they  are  related  by  the  definition 
of  a  unit  turbulent  Froudc  number,  that  is, 
-  1,  where  N  is  the  Brunt-VfiislilS  fre¬ 
quency.  taken  to  be  I0~^  s~'.  The  horizontal  turbu¬ 
lence  scale  in  the  boundary  layer  is  taken  to  be  the 
same  as  in  the  free  atmosphere,  as  discussed  below. 

Wind  profiles  were  obtained  directly  from  the  NGM 
fields  with  no  boundary-layer  correction  imposed.  The 
scheme  suggested  by  Draxler  ( 1990)  implies  a  linear 
velocity  profile  below  the  first  NGM  g^  point  (roughly 
ISOm)  under  stable  conditions.  This  was  found  to  pre¬ 
dict  unrealistically  low  translation  speeds  for  the  low- 
level  releases  in  ANATEX.  where  extensive  stable  con¬ 
ditions  were  predicted.  Much  better  predictions  were 
obtained  with  no  correction  and  by  simply  using  the 
lowest  NGM  grid-level  velocity  at  all  points  between 
that  level  and  the  surface.  Bas^  on  the  model  results 
described  later,  we  suggest  that  the  NGM  boundary 
layers  may  be  biased  toward  shallow  stable  layers,  es¬ 
pecially  in  January  and  February,  where  extensive  areas 
of  stable  conditions  were  predicted  for  continuous  pe¬ 
riods  of  many  days.  Some  calculations  were  therefore 
run  with  a  “standard”  boundary  layer,  which  was  uni¬ 
formly  imposed  over  the  entire  domain  and  exhibited 
the  diurnal  variation  illustrated  in  Fig.  2.  The  height 
variations  and  surface  heat  flux  were  specified  as  shown, 
and  the  local  time  used  for  the  variation  was  CST  over 
the  entire  domain.  The  numerical  values  used  in  the 
runs  were  hm,,  =  1000  m  for  the  January-February 
period,  and  hmn  =  1 500  m  for  March,  and  hmn  =  50  m 
for  both  periods.  The  surface  heat  flux  associated  with 
this  boundary  layer  is  zero  during  the  night  and  a  cosine 
variation  during  the  day.  with  a  maximum  value  of 
0.05®C  m  s‘‘ .  This  is  clearly  a  very  approximate  rep¬ 
resentation  and  cannot  represent  any  particular  day. 
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Fic.  2.  Standard  boundary-layer  diurnal  vanation. 


but  it  provides  a  reasonable  average  diurnal  variation 
of  the  boundary  layer. 

d.  Horizontal  wind  fluaucuions 

The  specification  of  dispersion  in  terms  of  velocity 
fluauation  statistics  avoids  the  use  of  completely  em- 
pirical  diSiisivities  but  requires  that  these  statistics  be 
provided  as  part  of  the  meteorolagical  description.  This 
can  be  view^  as  a  more  complete  specification  of  the 
flow  conditions,  defining  the  unresolved  or  random 
component  of  the  velocity  in  addition  to  the  mean 
wind.  Two  recent  pieces  of  research  have  allowed  us 
to  characterize  this  component  of  the  large-scale  flow 
in  a  nonarbitrary  way,  although  not  with  the  level  of 
understanding  that  we  have  for  small-scale,  three-di¬ 
mensional  turbulence  in  the  planetary  boundary  layer. 
First.  Nastrom  and  Gage  ( 1985)  analyzed  the  GASP 
(Global  Atmospheric  Sampling  Program)  wind  data, 
collected  by  specially  instrumented  commercial  aircraft 
over  a  long  period.  This  work  delineated  the  spearal 
behavior  of  the  wind  fluauations  in  the  upper  tropo¬ 
sphere  and  lower  stratosphere  as  a  funaion  of  latitude 
and  season.  The  data  shows  a  clear  spectral  break  at  a 
horizontal  wavelength  of  about  400  km.  with  a  k~^ 
spectrum  at  longer  wavelengths  and  for  shorter 
scales.  Second.  Gifford  ( 1988)  proposed  a  theoretical 
model  to  describe  the  spectrum,  arguing  that  the  large- 
scale  bduvior  is  consistent  with  the  predictions  of  two- 
dimensional  turbulence  theory,  while  the  behav¬ 
ior  indicates  three-dimensional  motions  in  that  part  of 
the  spectrum.  Gifford's  analysis  assumes  an  integral 
time  scale  of/~'  for  the  three-dimensional  turbulence, 
where  /  is  the  Coriolis  parameter,  using  dispersion 
measurements  as  supporting  evidence  ( Gifibrd  1 985 ). 
A  corresponding  turbulence  length  scale  is  easily  de¬ 
rived  from  the  time  and  velocity  scales. 

We  assume  that  the  upper-troposphere  data  is  rep¬ 
resentative  of  the  entire  lower  atmosphere,  althou^ 
there  is  probably  some  vertical  variation  in  practice. 
At  present,  the  data  do  not  exist  to  define  the  fluau¬ 


ation  profiles,  and  the  simplest  assumption  is  made. 
Furthermore,  the  GASP  data  represent  long-term  av¬ 
erages  and  cannot  provide  any  relationship  between 
the  local  meteorological  conditions  and  the  local 
smaller-scale  velocity  fluauations.  Use  of  the  GASP 
dau  is  therefore  equivalent  to  a  climatologically  av¬ 
eraged  spearum.  Using  Gifford’s  suggestions  and  the 
Gasp  profiles,  which  show  an  increase  in  fluauation 
variance  toward  the  poles,  the  horizonul  velodty  vari¬ 
ance  is  represented  as 


»'g/(2costf).  e<75'’ 

r  5/(2  00575“).  »s^75°. 


where  is  the  latitude  and  I  '  5  =  5.5  m*  s  .  The  value 
of  K5  is  suggested  by  Gifford  ( 1988 )  using  an  energy 
dissipation  estimate,  but  the  inverse  cosine  variation 
is  an  empirical  faaor  to  provide  a  match  with  the  lat¬ 
itudinal  variation  observed  in  the  GASP  data.  The  tur¬ 
bulence  length  scale  associated  with  these  dissipative 
fluauations  is  taken  to  be  Ar  =  (u’f  +  VTV'/f.  The 
singularity  in  /  at  the  equator  is  not  a  problem  in  the 
current  applications  since  the  southern  boundary  of 
the  domain  is  20°  N. 

The  velocity  variance  used  in  the  dispersion  calcu¬ 
lation  depends  on  the  resolution  of  the  input  wind  fields 
since  the  turbulence  input  represents  only  the  unre¬ 
solved  component  of  the  wind.  Using  the  filter  scale. 
A<;,  which  is  related  to  the  wind-field  resolution,  we 
define 


A//  *  minfAc,  Ar) 


consistent  with  a  spectral  behavior.  For  the  winds 
used  in  this  study,  we  postulate  Ac  >  Ar  since  Ar  30 
km  and  the  numerical  grid  length  is  about  1()0  km. 
The  model  calculations  will  therefore  use  the  full  spec¬ 
tral  energy,  u'f. 

This  essentially  specifies  all  the  horizontal  turbulence 
statistics  required  for  a  dispersion  prediction,  and  the 
results  of  the  calculations  using  these  statistics  will  be 
presented  in  section  4. 


e  Model  performance  measures 

It  is  well  known  that  trajeaory  errors  are  one  of  the 
major  causes  of  poor  conation  between  model  pre¬ 
dictions  and  observations.  An  error  in  plume  location 
will  result  in  a  negative  correlation,  that  is.  zero  pre¬ 
diction  for  samplers  impaaed  by  the  plume  and  high 
concentrations  predict^  where  nothing  is  observed. 
The  errors  in  trajeaory  prediction  are  largely  due  to 
wind  field  errors,  although  some  discrepancy  can  result 
from  incorrea  prediaion  of  the  vertical  distribution 
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of  polluunts.  The  errors  in  the  NGM  trajectory  pre¬ 
dictions  are  discussed  by  Draxler  ( 1991 )  and  in  the 
model  comparison  study  of  Cark  and  Cohn  ( 1990) 
and  will  not  be  present^  here.  Suffice  it  to  say  that 
the  NGM  errors  lie  within  the  range  of  the  other  tra¬ 
jectory  methods  reported  by  Clark  and  Cohn. 

We  will  concentrate  on  the  statistical  distribution  of 
the  concentration  observations,  uncorrelated  in  time 
and  space.  Some  localization  will  be  introduced  by 
consideration  of  two  separate  time  periods  of  obser¬ 
vations  and  also  by  partition  by  distance  from  the 
source.  Clearly,  as  the  sampled  set  is  made  more  spe¬ 
cific.  the  number  of  observations  are  reduced  and  the 
statistical  comparison  becomes  less  significant.  The 
relatively  large  sample  size  for  the  surface  observations 
allows  a  reasonable  partition  to  be  made,  and  cumu¬ 
lative  concentration  distributions  will  be  presented  for 
a  selection  of  samples.  Similar  comparison  will  be  made 
for  the  aircraft  data. 

4.  Results 

a.  Ground  samplers 

The  three-month  ANATEX  study  is  divided  into 
two  parts  by  a  missing  week  of  NGM  wind  fields,  firom 
23  February  to  1  March.  This  gives  a  seven-week  con¬ 
tinuous  period  from  S  January  to  21  February,  and  a 
four-week  period  in  March.  Since  the  model  needs  to 


be  reinitialized  in  March,  several  days  are  allowed  for 
the  tracers  to  become  esublished  on  the  sampler  net¬ 
work  and  comparisons  in  the  second  period  are  re¬ 
stricted  to  cover  9  March  to  28  March.  This  is  not  the 
same  division  as  used  by  Clark  and  Cohn  ( 1990),  who 
considered  two  six-week  periods,  but  is  the  most  con¬ 
venient  means  of  treating  the  missing  wind  data. 

Typical  maps  of  the  predicted  instantaneous  con¬ 
centrations  at  two  altitudes  are  shown  in  Figs.  3  and 
4.  These  maps,  which  show  tracer  patterns  from  many 
releases  of  PTCH  from  Glasgow  and  oPDCH  from  St. 
Cloud,  illustrate  the  different  characteristics  of  the  two 
periods.  The  January-February  winds  are  generally 
northwesterly  with  a  very  shallow  mixed  layer,  resulting 
in  very  little  tracer  being  found  at  1  km  in  Fig.  3.  The 
vertic^  wind  shear  is  also  evident  in  Fig.  3,  with  the 
tracer  at  1  km  displaced  to  the  east  of  the  surface-level 
tracer.  Figure  4  shows  a  typical  March  distribution, 
with  very  light  winds  and  dispersion  over  the  entire 
domain.  The  deeper  mixing  is  also  evident  with  almost 
uniform  concentrations  up  to  1  km. 

The  surface  sampler  network  records  24-h  avera^ 
concentration,  and  the  time-dependent  SCIPUFF  pre¬ 
diction  must  be  averaged  over  the  appropriate  period 
for  comparison  with  the  observations.  A  typical  time 
history  for  a  sampler  in  the  center  of  the  network  is 
shown  in  Fig.  S,  along  with  the  24-h  average  and  the 
observed  concentrations.  The  individual  releases  are 
often  evident  in  these  traces,  although  the  time  of  pas- 
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Fig.  3.  SCIPUFF  instantaneous  PTCH  and  oPDCH  surface  concentration  patterns 
at  two  altitudes  for  1800  UTC  29  January. 
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Fig.  6a.  Comparison  of  the  SCIPUFF  24-h  averase  PTCH  surtace  concentrauon  patter;  with  the  observed  patterns 

for  the  penod  28-30  Januarv-. 


sage  varies  considerably,  and  it  is  sometimes  difficult 
to  distinguish  between  releases.  There  are  periods  of 
good  correlation  between  observed  and  predicted  con¬ 
centration.  but  in  general  there  is  a  poor  prediction  of 
the  actual  value  on  a  day-by-day  basis.  This  is  consistent 
with  other  model  comparisons,  as  discussed  in  the  In¬ 
troduction.  and  is  largely  attributable  to  errors  in  flow 
trajeaory. 

In  addition  to  trajeaory  errors,  there  is  an  inherent 
uncertainty  due  to  unresolved  turbulence.  SCIPUFF 
contains  a  prediaion  of  the  turbulence-induced  uncer¬ 
tainty.  but  the  effea  of  24-h  time  averaging  (discussed 
in  appendix  B )  reduces  the  fluauation  intensity  to  near- 
negligible  proportions.  Therefore,  only  mean  concen¬ 
tration  is  considered  in  the  analysis  of  the  24-h  average 
dau.  The  turbulent  fluauation  component  will  be 
discussed  later  in  conneaion  with  the  short-term  air¬ 
craft  data. 


Figure  6  shows  two  3-day  sequences  of  the  24-h  av¬ 
erage  surface  tracer  concentrations  compared  with  the 
equivalent  maps  from  the  observed  values.  In  general, 
the  shape  and  location  of  the  prediaed  surface  con¬ 
centration  patterns  show  good  overall  agreement  with 
the  observations  for  both  periods.  This  is  an  indication 
of  the  accuracy  of  the  prescribed  wind  field  and.  also, 
the  lateral  diffusion  parameterization.  The  concentra¬ 
tion  magnitude  is  evidently  overprediaed  in  the  Jan¬ 
uary  period:  however,  the  quantitative  agreement  is 
much  bater  for  the  March  period.  This  discrepancy 
will  be  discussed  further  in  conneaion  with  the  statis¬ 
tical  analysis  below. 

The  trajeaory  errors  have  been  analyzed  by  Qark 
and  Cohn  ( 1990)  for  a  range  of  wind  models,  including 
the  NMC  NGM  winds,  so  no  such  analysis  is  performed 
here.  We  proceed  direaly  to  a  comparison  between  the 
cumulative  distributions  of  the  sampler  concentrations. 
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Fig.  6b.  Companson  of  the  SCIPUFF  24-h  avenge  PTCH  surface  concentntion  patterns  with  the  observed  patterns 

for  the  period  23-25  March. 


These  distributions  include  all  the  samplers  over  the 
whole  time  period,  with  the  exception  of  the  three 
samplers  closest  to  St.  Cloud  for  PMCH  and  oPDCH 
and  the  closest  sampler  to  Glasgow  for  PTCH.  These 
samplers  are  within  the  puff-initialization  region,  and 
therefore,  accurate  predictions  cannot  be  made  for 
them.  Figure  7  shows  the  cumulative  distribution  for 
the  January-February  period  and  indicates  an  over- 
p^ction  by  roughly  a  factor  of  2  for  all  three  tracers. 
Since  the  horizontal  spread  and  location  of  the  surface 
concentration  patterns  are  generally  well  predicted,  as 
illustrated  in  Fig.  6.  the  most  likely  cause  of  the  over- 
prediaion  is  incorrea  vertical  diffusion  of  the  tracer. 
An  examination  of  the  NGM  boundary-layer  predic¬ 
tion  shows  that  persistent  shallow,  stable  layers  occur 
frequently  in  the  January-February  period  over  much 
of  the  domain.  In  contrast,  the  March  prediaions  typ¬ 
ically  show  deeper  conveaive  layers  developing  almost 


every  day.  We.  therefore,  suggest  that  the  overpredic¬ 
tion  may  be  due  to  insufficient  vertical  mixing  using 
the  NGM  boundary-layer  predictions.  Figure  7  also 
shows  the  results  using  a  standard  1000-m-deep  mixed 
layer  every  day  over  the  whole  spatial  domain  for  the 
January-February  period.  The  dispersion  prediction  is 
obvioudy  sensitive  to  the  boundary-layer  estimate  since 
this  determines  the  vertical  dilution  of  the  tracer:  a 
1000-m  layer  apparently  gives  a  much  better  distri¬ 
bution.  The  results  can  be  examined  in  more  detail  by 
subdividing  the  samplers  into  arcs  depending  on  their 
distance  from  the  sources.  The  distributions  for  PTCH 
on  three  separate  arcs,  shown  in  Fig.  8.  provide  evi¬ 
dence  that  the  model  is  predicting  the  correct  evolution 
of  the  tracer  with  distance. 

Similar  results  are  available  for  the  three-week  March 
period,  although  use  of  the  NGM  boundary  layer  gives 
a  better  prediction  here,  as  shown  in  Fig.  9.  The  results 
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Fic.  7.  Predicted  and  obseryed  ground  sample  cumulative  distri¬ 
bution  functions  for  the  period  5  January-21  Febniary  < - ob¬ 
served;  predicted;  NGM  boundary  layer - .  standi  lOOO-m 

boundary  layer . ). 


for  a  sundard  I  S00>m  layer  are  also  presented  in  Fig. 
9  and  show  an  improvement  The  distributions  on  the 
arcs  show  a  good  prediction  of  the  downwind  evolution 
for  this  three-week  period,  but  it  is  similar  to  Fig.  8 
and  is  therefore  not  shown. 

b.  Aircraft  data 

The  data  from  the  aircraft  samplers  provide  short¬ 
term  spatial  averages  over  a  pathlength  typically  on  the 
order  of  20  km.  The  corresponding  model  predictions 
are  computed  as  the  average  of  the  mean  concentration 
at  the  flight  path  endpoints  ( effectively  assuming  a  lin¬ 
ear  variation  with  path  length ).  Reductions  in  the  pre¬ 
dicted  concentration  variances  due  to  spatial  averaging 
are  estimated  by  a  procedure,  describe  in  appendix 
B.  based  on  an  assumed  spatial  correlation  function 
with  an  integral  scale  related  to  the  internal  fluctuation 
scale  Ih- 


Turbulent  fluauations  are  much  more  significant  in 
these  shorter-term  averages  than  in  the  24-h  surface 
averages,  and  any  comparison  with  the  dau  requires 
that  account  of  the  resulting  uncertainty  be  made.  Fol¬ 
lowing  the  philosophy  of  Lewellen  and  Sykes  ( 1989). 
the  aircraft  dau  are  therefore  treated  as  r^izations  of 
a  random  variable  theoretically  defined  by  appropriate 
probability  distribution  funaions.  Our  modeling  as¬ 
sumes  a  form  for  the  distribution  funaion  that  can  be 
defined  using  only  the  predicted  mean  and  variance. 
As  in  our  shorter-range  results,  we  assume  a  clipped 
normal  distribution,  since  it  gives  better  agreement  than 
the  lognormal  distribution  suggested  by  Csanady 
( 1973).  The  clipped  normal  distribution  resulu  from 
replacing  any  unrealizable  negative  tail  of  a  normal 
distribution  by  a  delu  function  at  zero  concentration 
whose  strength  is  defined  by  integrating  over  the  neg¬ 
ative  range.  Details  can  be  found  in  Lewellen  and  Sykes 
(1986). 


(Concentration  (flA) 


Fig.  8.  Prediaed  and  observed  ground  sample  cumulative  distri¬ 
bution  function  for  PTCH  along  three  arcs  for  the  period  5  January- 
2 1  February.  Predictions  were  made  using  a  standard  boundary-layer 
depth  of  1000  m  ( - observed: - predicted). 
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Concentration  ((LA.) 


Fig.  9.  Predicted  and  otMcrved  ground  jample  cumulative  distri¬ 
bution  functions  for  the  period  9-28  March  ( - observed;  predicted; 

NGM  boundary  layer - ,  standard  ISOO-m  bound^  layer 

. ). 


The  cumulative  distribution  of  aircraft  sampler  con¬ 
centrations  of  PTCH  from  Glasgow  is  shown  in  Fig. 
10.  The  data  is  restricted  to  samples  further  than  SO 
km  from  the  source  to  avoid  the  puff-initialization  re¬ 
gion.  This  is  a  minor  restriction  since  most  of  the  sam¬ 
ples  were  from  greater  disunces.  The  figure  also  shows 
the  distributions  of  the  predicted  mean  concentration 
c,  the  expected  distribution  with  90%  confidence  limits, 
and  the  observations.  The  expected  distribution  is  ob¬ 
tained  by  making  a  large  number  of  random  samples 
from  each  of  the  predicted  distributions.  Specifically, 
for  each  aircraft  sample,  the  model-prediaed  mean  and 
variance  ( of  the  average  along  the  flight  path )  are  used, 
along  with  the  assumed  clii^xd-normal  shape,  to  define 
a  probability  distribution.  A  set  of  pseudo-observations 
are  created  by  taking  random  samples  from  each  dis¬ 
tribution  and  constructing  the  resulting  cumulative 
distribution  function  (cdO-  After  taking  a  large  number 
of  random  samples  and  constructing  a  large  number 


of  cdf ’s  from  these  pseudo-observations,  the  expected 
distribution  is  obtained,  that  is,  the  average  c^  and 
confidence  bounds.  (The  90%  confidence  bounds 
shown  in  Fig.  10  indicate  that  90%  of  the  cdf’s  con¬ 
structed  from  individual  seu  of  the  pseudo-observa¬ 
tions  fall  between  the  two  appropriate  curves.) 

The  contrast  between  the  c  distribution  and  the  ex- 
peaed  distribution  is  worth  noting.  Here  c  is  the  pre- 
diaed  ensemble-average  plume  concentration,  taldng 
no  account  of  fluctuations,  and  therefore  tends  to  be 
spread  over  a  larger  region  than  that  occupied  by  an 
actual  plume  realization  (observed  or  “predicted”). 
Thus,  the  mean  distribution  contains  many  small  val¬ 
ues  but  few  nonzero  predictions,  while  the  expected 
distribution  contains  a  small  number  of  large  values 
but  many  zero  predictions,  that  is,  high  intermittency. 

Figure  10  illustrates  some  features  that  are  typical 
of  all  the  predicted  tracer  distributions.  The  distribution 
of  c  conuins  too  few  small  values  ( <  1  fl  F' )  and  too 
few  large  values  for  the  reasons  just  discussed,  while 
the  expected  distribution  has  too  many  small  values 
compared  with  the  observations.  The  discrepancy  be¬ 
tween  the  observed  and  expected  distribution  is  attrib¬ 
uted  to  the  nonrandom-sampling  procedure  used  in 
the  aircraft  sampling.  Flights  were  made  on  the  basis 
of  the  most  recent  wind  forecasts,  with  the  objective 
of  intersecting  the  tracer  plume.  Furthermore,  on-board 
real-time  sampling  was  used  to  locate  the  tracer  and 
to  determine  the  flight  pattern  during  the  mission  to 
provide  several  transects  of  the  plume.  This  procedure 
can  be  expected  to  produce  fewer  “zero”  samples  than 
a  completely  random  sampling,  and  some  allowance 
for  this  effect  needs  to  be  made  in  the  comparison. 

This  conditional  sampling  of  the  tracer  cannot  be 
modeled  easily,  and  we  therefore  arbitrarily  match  the 
observed  and  predicted  distributions  at  the  1  fl  1~'  level. 


Concentration  (fUL) 

Fig.  10.  Predicted  and  observed  aircraft  sample  cumulative  dis¬ 
tribution  functions  for  the  period  9-28  March.  Predictions  were  made 

using  a  standard  boundary-layer  depth  of  1 500  m  ( - observed; 

predicted:  expected - .  mean - .  90%  confidence  limits 

- >. 
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Fig.  II.  Predicted  and  observed  aircnfi  sample  cumulative  disnibution  functions  I - observed: - predicted: - 90%  confidence 

limits).  PredictMiis  were  made  using  a  standard  boundary  layer,  (a)  Period  S  January-21  February  with  lOOO-m  boundary  layer,  (b)  Period 
9-28  March  with  I  SOO-m  boundary  layer. 


This  level  was  chosen  as  being  effectively  zero  concen¬ 
tration  so  that  only  the  distribution  of  nonzero  con¬ 
centrations  can  be  compared,  that  is.  concentrations 
greater  than  I  fl  r' .  The  small  concentrations  are  ig¬ 
nored  since  the  fraction  of  “zero"  aircraft  observations 
depends  on  the  plume-locating  skill  of  the  flight  ob¬ 
servers. 

The  adjusted  results  for  the  two  periods  using  the 
standard  boundary-layer  e^mate  are  shown  in  Fig. 
1 1,  which  includes  the  expected  distribution  together 
with  90%  confidence  bounds.  The  comparison  between 
the  predicted  distributions  and  the  airciaft  cAiservations 
is  encouraging,  although  we  note  again  the  forced 
agreement  at  1  fl  1 ' ' .  The  information  provided  by  the 
variance  allows  us  to  give  a  much  better  estimate  of 
the  expected  short-term  concentrations  and  also  some 
measure  of  the  confidence  bounds  on  the  prediction. 
The  clipped  normal  function  provides  a  reasonable  de¬ 
scription  of  the  shape  of  the  cumulative  distribution. 


especially  for  the  March  period  where  the  observations 
lie  almost  completely  within  the  90%  confidence  bands. 

It  should  be  poinuxl  out  that  the  sampling  procedure 
just  described  ignores  correlation  betwera  observations 
since  each  random  choice  is  made  independently.  Since 
the  pathlengths  do  not  generally  provide  more  than  a 
few  observations  within  the  tracer  cloud  we  do  not 
have  the  situation  of  strong  correlation  between  adja¬ 
cent  samplers,  and  this  is  not  a  severe  limitation.  The 
effect  of  such  correlation,  in  general,  is  to  reduce  the 
number  of  degrees  of  freedom  in  the  problem  and 
hence  increase  the  width  of  the  confidence  band. 

c.  Ensemble  predictions 

The  partitioning  of  the  velocity  fluauations  into  re¬ 
solved  and  turbulent  components  depends  on  the  res¬ 
olution  of  the  measured  wind  field.  As  the  detail  in  the 
wind  field  input  is  reduced,  the  uncertainty  in  the  dis- 
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person  prediction  is  increased,  along  with  the  turbulent 
energy  component.  The  applicability  of  our  statistical 
prediction  can  be  investigated  by  reducing  the  reso¬ 
lution  of  the  wind  input  and  testing  the  prediaions 
based  on  the  lower  information  content.  An  extreme 
test  is  provided  by  using  an  average  of  the  wind  field 
over  the  entire  spatial  domain  and  over  the  two  time 
periods  of  interest.  The  resulting  wind  vanes  only  with 
altitude  and  represents  a  long-term  average  for  the 
eastern  United  States. 

The  mean  and  rms  fluauation  component  profiles 
for  the  two  periods  of  NGM  wind  data  are  shown  in 
Fig.  12.  The  mean  wind  is  generally  westerly  and  in¬ 
creases  with  height:  the  profiles  are  very  similar  except 
close  to  the  surface  where  the  March  data  shows  near- 
stagnant  conditions  in  the  average.  The  rms  fluauation 
profiles  are  also  very  similar,  with  a  value  of  6-7  m  s  ~ ' 
in  the  lowest  few  kilometeis  and  increasing  up  to 
around  IS  m  s~'  at  r  =  8  km. 

Idealized  profiles  were  specified  for  the  model  runs. 
Only  the  lowest  few  kilometers  are  relevant  since  the 
tracer  is  generally  contained  below  3  km.  Linear  ve¬ 
locity  profiles  were  fit  to  the  mean  velocity  and  a  con¬ 
stant  value  of  SO  m^  s~*  was  used  for  the  horizontal 
velocity  variances.  A  length  scale  for  the  horizontal 
velocity  fluauations  also  needed  to  be  specified,  and 
it  was  found  that  Am  ~  2000  km  gave  a  reasonably 
consistent  description.  First  this  scale  is  typical  of  the 
large-scale  eddies  in  the  NGM  field,  and  second  the 
assumption  of  the  k~^'^  spectrum  gives  smaller-scale 
energies  close  to  those  assumed  in  the  previous  cal¬ 
culations.  Thus,  if  the  typical  scale  of  Ar  *  30  km  from 
section  3d  is  used,  the  appropriate  energy  at  that  scale 
in  the  ensemble  flow  would  be  30(30/2000)^^^  as  3 
m^  s~^.  It  should  be  noted  that  the  energy  in  the  larger 
scales  is  probably  not  proportional  to  but  we 
believe  that  the  most  important  feature  here  is  the  rep- 
resenution  of  the  bulk  of  the  spearum.  Over  the  trans¬ 
port  times  of  this  study,  the  concentration  fluauation 
scales  do  not  grow  lai^r  than  2000  km  and  the  rep¬ 
resentation  of  the  long-wave  end  of  the  spearum  is  not 
critical.  The  specification  is  best  regarded  as  investi¬ 
gative  and  we  hope  to  show  that  further  study  is  war¬ 
ranted. 

The  dispersion  calculation  with  the  ensemble  wind 
field  is  very  simple,  since  the  wind  is  steady  and  ho¬ 
mogeneous.  The  standard  diumally  varying  boundary 
layers  are  used,  which  are  representative  of  an  ensemble 
average,  for  the  calculation  with  a  daytime  maximum 
depth  of  1000  m  for  January-February  and  1300  m 
for  March.  The  calculations  are  computationally  fast, 
using  only  a  few  puffs  that  spread  rapidly  with  the  high 
turbulence  levels  and  light  wind. 

We  can  compare  only  the  cumulative  distribution 
of  the  surface  samples,  since  the  prediction  uses  no 
information  about  particular  days  in  the  period.  The 
mean  concentration  prediaion  is  a  very  slow-moving, 
rapidly  spreading  circular  Gaussian  shape  with  the  day¬ 
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Fig.  12.  Ensemble  velocity  ptobics  ( - u: - v). 


March  4  06:00  -  March  29  24:00 


to-day  variability  represented  statistically  by  the  vari¬ 
ance  inediction.  Expected  cumulative  di^lwtions  are 
generated  in  the  same  way  as  the  aircraft  predictions: 
by  selecting  randomly  from  a  probability  distribution 
at  each  sampler  location.  The  variance  was  reduced  by 
the  appropriate  averaging  over  the  24-h  period  of  the 
observation,  as  described  in  appendix  B.  and  the  log¬ 
normal  probability  distribution  function  was  found  to 
give  a  much  better  fit  to  the  data.  The  intermittency 
of  the  clifiped-normal  distribution  gave  poor  agreement 
with  the  observations,  perhaps  as  a  result  of  the  long 
time  averapng  that  reduces  the  probability  of  a  zero 
concentration. 

The  results  of  the  simplified  ensemble-average  pre¬ 
diction  for  the  two  time  periods  are  shown  in  Figs. 
i3a.b.  The  prediaion  is  remarkably  good  considering 
the  limited  sutistical  information  used  in  the  calcu¬ 
lation.  The  use  of  the  concentration  variance  prediction 
and  the  assumed  lognormal  probability  distribution 
allows  an  estimate  of  the  24-h  ave  'zge  concentrations 
from  the  much  longer  averaging  times  in  the  wind  sta¬ 
tistics.  The  March  statistics  are  very  well  prediaed. 
except  for  the  lowest  concentrations  of  PTCH.  It  is  not 
clear  that  the  distribution  below  1  fl  r'  is  significant 
since  the  background  noise  levels  are  of  this  magnitude 
(see  Table  1 ).  but  the  results  down  to  0.1  fl  T'  are 
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hici .  13.  Predicted  and  observed  ground  sample  cumulative  distribution  functions.  Prediaions  were  made  using  the  ensemble  meteoroiofy 
( - observed: - predicted; . 90%  confidence  limits),  (a)  Period  5  Januarv-2!  February,  (b)  Period  9-28  March. 


shown  for  completeness.  The  most  significant  error  is 
in  the  PTCH  distribution  in  the  January-February  pe¬ 
riod.  which  shows  an  underprediction  of  the  values 
above  1  fl  r' .  The  other  tracers  are  better  predicted, 
however,  and  generally  show  much  closer  resemblance 
to  the  March  period  distributions.  The  model  results 
are  also  very  similar  for  the  two  periods,  since  the  wind 
statistics  do  not  differ  significantly. 

Some  skill  is  demonstrated  in  the  prediction  of  the 
variation  with  distance,  as  shown  in  Fig.  14.  Here  the 
cumulative  distributions  of  PTCH  samples  for  March 
are  shown  along  various  sampler  areas.  We  choose  the 
March  period  because  of  the  reasonably  good  overall 
prediction,  and  we  have  combined  adjacent  areas  to 
increase  the  sample  size  and  reduce  the  confidence 
bounds.  The  range  of  the  confidence  bounds  is  still 
large,  but  the  model  does  show  the  general  decrease 
with  disunce  from  the  source. 


5.  Summary  and  cooclnsioas 

A  Lagrangian  puff  model  for  long-range  atmospheric 
dispersion.  SCIPUFF.  has  been  described.  The  model 
uses  second-order  turbulence  closure  methods  to  rep¬ 
resent  the  diffusion  processes  and  has  been  develop^ 
from  shorter-range  applications  where  boundary-layer 
turbulence  dominates.  The  mesoscale  velocity  fluctu¬ 
ation  statistics  for  longer-range  calculations  are  based 
on  recent  descriptions  of  the  energy  spectrum  from 
measurements  and  theoretical  considerations.  The  La¬ 
grangian  framework  provides  an  efficient  description 
of  the  local  concentration  field  and  is  coupled  with 
splitting  and  merging  algorithms  to  maintain  an  eco¬ 
nomical  representation  with  appropriate  resolution  of 
the  flow  features.  The  internal  puff  concentration  dis¬ 
tribution  is  a  fully  general  Gaussian  shape  so  that  shear 
distortion  effects  can  be  computed  accurately. 
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Fig.  14.  Predicud  and  otMerved  ground  sample  cumulative  dis¬ 
tribution  function  for  oPDCH  along  three  pairs  of  arcs  for  the  petiod 
9-28  March.  Predictions  were  made  using  the  ensemble  meteorolagy 
( - observed;  -  -  -  predicted: . 90%  confidence  limits). 


The  turbulence  closure  includes  a  prediction  for  the 
concentration  fluauation  variance,  a  measure  of  the 
random  variability  in  an  observable  value.  This  can  be 
used  to  provide  a  measure  of  the  statistical  significance 
in  model  evaluation.  The  variance  prediaion  relates 
the  uncertainty  in  the  input  wind  field,  that  is.  the  tur¬ 
bulence  velocity  fluctuations,  to  the  uncertainty  in  the 
predicted  concentration  field.  Turbulent  fluauations 
are  generally  considered  to  be  relatively  small  in  scale, 
but  the  spectrum  of  atmospheric  motions  is  continuous 
and  the  definition  of  the  turbulent  wind  component 
depends  on  the  detail  of  the  available  resolved  wind 
field.  One  could  also  consider  input  wind  inaccuracy 
as  a  stochastic  component  with  a  measurable  spectrum, 
but  our  present  understanding  of  the  large-scale  wind 
fluauations  is  not  sufficiently  well  developed  for  a 
quantitative  assessment  of  the  various  contributing 
components. 


ET  AL. 

The  use  of  explicit  velocity  fluctuation  statistics,  and 
consideration  of  the  spectral  distribuuon  of  the  fluc¬ 
tuation  energy,  allows  a  rational  approach  to  the  pre¬ 
diction  of  observable  concentrations  and  the  effects  of 
sampling  times,  in  contrast  to  the  use  of  an  empirical 
dispersion  parameter  based  on  averaging  time,  the  tur¬ 
bulence  closure  framework  provides  an  ensemble  de¬ 
scription  of  the  concentration  field  and  a  method  for 
predicting  the  statistics  of  any  sampling  procedure.  The 
prediaion  requires  assumptions  for  the  probability 
distribution  funaions  and  correlation  funaions.  but 
these  can  be  made  on  a  rational  basis  with  the  potential 
for  direa  measurement  of  these  quantities  at  some  ume 
in  the  future. 

The  SCIPUFF  model  has  been  exercised  on  the  AN- 
ATEX  dispersion  dau  using  the  wind  fields  from  the 
NMC  Nested  Grid  Model.  The  experiment  recorded 
three  months  of  24-h  average  tracer  concentrations  on 
a  sampla  nawork  ranging  3(XX)  km  from  the  source. 
The  impaa  from  multiple  short-duration  releases  from 
two  separate  sources  was  measured  at  77  ground-based 
samplers,  providing  a  large  database  for  evaluation 
purposes.  In  addition,  aircraft  samples  of  short  duration 
were  obtained  within  a  few  hundred  kilomaers  of  the 
sources.  The  model  predictions  and  observed  concen¬ 
trations  cannot  be  compared  on  a  paired  time  and  lo¬ 
cation  basis  due  in  part  to  the  inaccuracy  of  the  wind 
field,  but  the  model  predictions  of  the  sutistical  dis¬ 
tributions  of  the  observed  concentrations  are  generally 
good.  There  is  a  tendency  for  the  model  to  overptedia 
the  surface  impaa  in  the  January-February  period, 
and  this  is  ascribed  to  the  frequent  shallow,  stable 
boundary-layer  predictions  of  the  NGM.  Use  of  a  sun- 
dard  lOGO-m-d^  daytime  layer  improved  the  predic¬ 
tion  sgnificantly.  A  reasonably  good  prediction  of  the 
aircraft  statistical  distribution  was  ^so  obtained  by 
matching  the  observed  ftaaion  of  zero  concentrations 
to  account  for  the  conditional  nature  of  the  aircraft 
sampling  and  assuming  a  clipped  normal  probability 
distribution  funaion. 

As  a  further  test  of  the  sutistical  prediction,  calcu¬ 
lations  were  made  with  the  average  wind  sutistics  for 
the  entire  NGM  domain  over  the  ANATEX  periods. 
This  modifies  the  turbulent  wind  component  to  include 
all  the  resolved  variations  in  the  NGM  wind  fields.  The 
predicted  concentration  variance  is  much  larga  in  this 
case,  and  a  reasonable  prediction  of  the  observed  dis¬ 
tribution  was  found  with  a  lognormal  probability  dis¬ 
tribution  for  the  24-h  samples. 

In  summary,  a  sutistical  dispersion  model  has  been 
compared  with  the  ANATEX  observations  using  two 
different  levels  of  resolution  in  the  wind  field.  The  re¬ 
sults  are  very  encouraging  but  would  benefit  greatly 
from  further  study  of  the  sutistical  nature  of  the  at¬ 
mospheric  flow  fields.  A  number  of  speculative  as¬ 
sumptions  were  necessary  in  the  application  of  the 
model  to  the  ANATEX  cases,  since  reliable  informa- 


108 


946 


JOURNAL  OF  APPLIED  METEOROLOGY 


Volume  32 


tion  on  spectral  intensities,  time  and  space  correlation 
scales,  and  probability  distribution  functions  is  not 
currently  available.  The  utility  of  a  statistical  prediction 
has  been  demonstrated  in  the  calculation  of  confidence 
bounds  in  the  comparison  with  the  short-term  samples 
from  aircraft  and  in  the  ensemble  wind-held  prediction. 
The  inherent  uncertainty  associated  with  any  predic¬ 
tion  of  atmospheric  transport  over  a  long  range  is  likely 
to  always  be  signihcant,  and  research  efforts  can  be 
devoted  to  describing  or  reducing  it.  The  results  pre¬ 
sented  above  suggest  that  a  quantitative  description  of 
the  uncertainty  is  within  reach. 
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APPENDIX  A 
Source  Specification 

In  principle,  the  plume  emanating  from  a  source 
with  a  mass  flow  rate  F  can  be  modeled  by  creating 
new  puffs  at  each  time  step  defined  by  a  mass  Q 
=  FA/,  a  centroid  defined  by  the  source  location  and 
second  moments  related  to  the  source  size.  The  sub¬ 
sequent  plume  evolution  would  then  simply  be  com¬ 
puted  using  SCIPUFF.  If  the  time  step  is  long  enough 
for  the  plume  to  travel  more  than  the  source  diameter, 
however,  then  the  spatial  evolution  of  the  plume  during 
the  time  step  must  be  considered.  Since  the  dispersion 
time  step  is  IS  min  and  the  source  size  is  less  than  I 
m.  plume  evolution  during  a  time  step  is  certainly  a 
consideration.  We  therefore  require  a  plume  model  to 
initialize  puffs  at  the  source  for  use  in  SCIPUFF.  A 
simple  m^el  for  the  plume  spread  might  be  adequate 
for  defining  mass,  location,  and  second  moments,  but 
the  rapid  evolution  of  the  fluctuation  variance  c'~ 
requires  more  sophisticated  treatment.  Therefore,  a 
plume  version  of  the  full  puff  equations  of  section  3 
are  solved:  that  is,  we  march  downstream  from  the 
source  with  a  very  small  time  step  for  a  time  /,  com¬ 
puting  dispersion  in  the  vertical  and  lateral  directions. 
The  source  time  i,  is  chosen  so  that  the  puff  time  scales 
are  large  enough  to  be  resolved  by  the  SCIPUFF  cal¬ 
culation  time  step.  A/;  /,  =  60  min  in  the  ANATEX 
calculations. 

At  each  time  step,  a  puff  is  initialized  at  the  wind- 
translated  location  corresponding  to  a  time  /,  after  re¬ 
lease.  with  a  streamwise  spread  of  UAl.  where  U  is  the 
local  wind  speed.  All  other  moments  are  obtained  from 
the  plume  calculation.  Thus,  for  example,  the  plume 
initialization  solves 

^  =  U(x„/)  (Al) 

at 


^  Fy,  =  2v?  -  AqHfy,,  (A3) 

at 

where  x,  is  the  source  location.  U(Xj,  /)  is  the  velocity 
at  the  source,  and  v','  is  the  velocity  variance.  Equations 
( Al )-( A3)  are  integrated  out  to  t  =  then  the  puff 
spatial  moments  are  rotated  from  the  local  frame  de¬ 
termined  by  the  plume  into  the  SCIPUFF  coordinate 
system.  The  two  horizontal  diffusivities.  Fi  u  and  F22a, 
are  initialized  with  the  value  of  F,,  from  the  plume 
calculation. 

The  calculation  of  the  mean-square  concentration. 
S„,  proceeds  similarly,  except  that  there  is  only  a  self¬ 
interaction  term  in  the  plume  equivalent  of  ( 1 8 )  and 
the  overlap  integral  is  two-dimensional.  These  plume 
equations  are  essentially  those  described  by  Sykes  et 
al.  ( 1986 ).  except  for  one  subtle  point  associate  with 
the  anisotropic  treatment  of  the  turbulence.  The  special 
treatment  of  the  vertical  dispersion  and  the  extra  dis¬ 
sipation  term.  (23).  are  retained  in  the  plume  calcu¬ 
lation.  This  implies  that  the  variance  production  is  only 
a  result  of  the  lateral  dispersion,  and  this  one-dimen¬ 
sional  situation  requires  that  the  coefficient  of  2bs  in 
(22)  he  reduced  to  bs  in  order  to  be  consistent  with 
the  Gifford  meandering  plume  result  (Gifford  1959). 


APPENDIX  B 


Effect  of  Averaging  on  the  Fluctuations 

The  field  observations  of  tracer  concentrations  in¬ 
volve  averaging,  either  purely  temporal  in  the  case  of 
the  24-h  ground  samples  or  mostly  spatial  for  the  air¬ 
craft  data.  SCIPUFF  predicts  the  instantaneous  fluc¬ 
tuations  variance,  which  is  reduced  by  the  averaging 
process,  and  the  reduction  must  be  accounted  for  before 
comparing  with  the  observations. 

If  the  average  concentration  is  defined  as 

I  rie*T/2 

r  =  -  c{t)dt,  (Bl) 

T  Jio-T/2 


then 

1  rio*TI2 

f  =  i  c(t)dt  (B2) 

T  Jit-T/2 


and 


r-2 


■pfj  c'(t)c'(t')dtdt’ 


where 


(B3) 


p(/;/')  = 


c'(t)c’(t') 

P(0 
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If  c(  / )  is  a  sutionary  process,  then  ( B3 )  reduces  to  the 
standard  result  ofTennekes  and  Lumley  (  1977 ):  that 
is. 


The  estimate  of  the  variance  in  an  averaged  obser¬ 
vation  requires  an  estimate  of  the  correlation  function. 
For  the  sutionary  case,  an  exponential  form  is  as¬ 
sumed,  so  that 


=  2c^ 


(/3-  1  +<•-") 


where  0  *  T/Tcind  pit')  = 

For  the  general,  nonsutionary  case,  a  crude  approx¬ 
imation  for  ( B3 )  is  used,  eifeaively  assuming 


Pin 


1.  U'l<n 

0.  |/'|>7-,. 


Then  ( B3 )  becomes 


r«2 


c^(r)min(/.  TcU)]dt. 


(B5) 


All  of  the  above  discussion  applies  immediately  to 
the  case  of  spatial  averaging  if  r  is  regarded  as  disunce 
along  the  path  of  integration  and  7,  is  replaced  by 
the  integral  scale  of  the  spatial  correlation. 

The  integral  scales  are  specified  in  terms  of  the  pre¬ 
dicted  fluctuation  scale  In  using  the  relations 


and 


r,= 


Ic  =  2/« 

(B6) 

0.75/« 

msK{L’.q„)  ■ 

(B7) 

The  time  scale.  (B7).  is  based  on  short-range  esti¬ 
mates  from  large-eddy  simulations  (Sykes  and  Henn 
1 992 ).  while  the  spatial  scale  is  determined  from  ideal¬ 
ized  calculations  for  randomly  displaced  Gaussian 
puffs.  Both  estimates  should  be  regarded  as  very  spec¬ 
ulative  at  this  stage  and  need  more  direa  measurements 
of  the  fluctuations  on  the  larger  scale. 
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